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Abstract 


A  three-dimensional  prism  acoustic  finite  element  was 
developed  to  model  cavity  resonances  of  three  axisymmetric 
combustion  chamber-cylinder  configurations.  The  lowest 
cavity  resonances  were  examined  in  detail  as  a  function  of 
cylinder  depth. 

Experimental  models  were  constructed  to  compare 
experimental  results  to  the  finite  element  results  for  both 
cavity  resonances  and  mode  shapes  for  each  configuration.  In 
addition,  analytical  theory  for  a  cylindrical  enclosure  was 
examined  in  some  detail.  Agreement  between  experimental, 
finite  element  and  analytical  results  was  good  for  the 
lowest  mode. 
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1.  INTRODUCTION 


1.1  Engine  Knock  Related  to  Acoustics 

Knock  in  the  internal  combustion  engine  is  a  phenomenon 
which  has  received  considerable  attention  in  the  past. 

Engine  knock  is  described  as  gas  vibrations  caused  by  a 
rapid  pressure  rise  in  the  combustion  chamber  of  the  engine 
(1,2,3).  An  audible  noise  or  "ping"  emanates  from  the  engine 
as  a  result  of  the  vibrations  of  the  gas-cylinder  system  at 
frequencies  of  the  order  of  5,000  Hz.  A  detailed  discussion 
of  engine  knock  is  given  by  Obert  (2).  He  describes  several 
of  the  undesirable  characteristics  of  engine  knock,  such  as 
engine  component  vibration,  scrubbing  action  on  the  chamber 
walls  and  localized  overheating,  which  suggest  that  knock 
should  be  avoided  in  the  combustion  engine.  Knock  prevention 
has  generated  a  need  to  determine  some  of  the  theoretical 
characteristics  of  the  knock  phenomenon. 

C.S.  Draper  (1)  proposed  a  theory  that  engine  knock  was 
related  to  the  lowest  acoustical  cavity  resonance  of  the 
combustion  chamber.  Draper  considered  an  engine  which 
consisted  of  a  circular  cylinder  with  flat  ends  at  right 
angles  to  the  cylinder  axis.  His  results  indicated  that  the 
knocking  frequency  of  the  engine  was  the  fundamental  cavity 
resonance  of  the  combustion  chamber.  This  was  an  important 
conclusion  as  it  provided  a  means  to  predict  the  frequency 
of  engine  knock  by  calculating  the  acoustical  cavity 


1 


. 


. 


2 


resonances  of  the  combustion  chamber-cylinder  configuration. 
This  information  is  required  for  design  of  suitable  knock 
detection  systems  (4). 

The  results  presented  by  Draper  (1)  were  based  on  the 
analytical  solution  of  the  wave  equation  for  a  cylindrical 
enclosure  as  given  by  Morse  (5).  Many  combustion  chambers 
are  not  simple  cylindrical  enclosures.  The  Internal 
Combustion  Engine  (6),  shows  a  selection  of  various 
combustion  chambers.  Often,  the  diameter  of  the  combustion 
chamber  differs  from  that  of  the  engine  cylinder.  As  knock 
occurs  near  top  dead  center  (2),  see  Figure  1.1,  and 
continues  for  a  portion  of  the  stroke,  the  combined  geometry 
of  the  combustion  chamber  and  cylinder  at  various  piston 
positions  (i.e.  cylinder  depths)  must  be  studied  to 
completely  describe  the  engine  knock  phenomenon.  A  study  on 
knock  in  spark  ignition  engines  (3),  mentions  the  dependence 
of  the  knocking  frequency  on  piston  position  shown  by 
earlier  investigators  (7). 

The  study  of  combustion  chamber-cylinder 
configurations,  such  as  those  shown  in  reference  (6)  is 
extremely  difficult  using  the  method  developed  by  Draper. 

The  difficulty  stems  from  the  irregular  boundary  of  the 
combustion  chamber-cylinder  configuration.  The  solution 
provided  by  Draper  was  for  an  axisymmetric  cylinder  solved 
analytically  using  the  wave  equation  in  cylindrical 
co-ordinates  with  the  given  boundary  conditions.  There  are 
inherent  difficulties  in  describing  a  combustion  chamber 
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-cylinder  configuration  if  the  diameter  of  the  combustion 
chamber  differs  from  that  of  the  cylinder.  A  technique  which 
lends  itself  to  solving  these  types  of  problems  is  the  use 
of  acoustic  finite  elements,  as  described  by  A.  Craggs  (8) 
to  obtain  an  approximate  numerical  solution. 

From  the  standpoint  of  cavity  resonances,  diesel  and 
spark  ignition  engines  are  treated  in  the  same  manner,  by 
considering  the  geometry  of  the  combustion  chamber  and 
cylinder  together.  Because  the  cavity  resonances  are  related 
to  the  geometry  considered,  this  analysis  is  concerned 
mainly  with  the  geometry  of  the  combustion  chamber-cylinder 
configurations  rather  than  the  engine  operating  conditions. 
Consequently,  the  study  considers  air  in  enclosures  at  room 
temperature . 

1.2  Acoustic  Finite  Element  Modelling 

Acoustic  finite  elements  were  presented  in  their 
simplest  form  by  G.M.  Gladwell  (9).  A  one-dimensional  column 
of  air  with  rigid  ends  was  examined  to  evaluate  the  natural 
frequencies  using  the  finite  element  method.  The  finite 
element  approximations  of  the  natural  frequencies  for  the 
one-dimensional  problem  were  in  good  agreement  with  the 
acoustic  theory.  This  spurred  further  examination  into  the 
usage  of  acoustic  finite  element  modelling. 

The  acoustic  finite  element  method  was  extended  to 
three-dimensional  problems  by  A.  Craggs  (8).  Different 
acoustic  finite  elements  were  used  to  determine  the  natural 


5 

frequencies  and  mode  shapes  of  complex  shaped  enclosures. 

Good  agreement  was  shown  between  the  finite  element 
predictions  and  experimental  results. 

Because  of  the  relationship  of  engine  knock  to 
acoustical  cavity  resonances  (1),  acoustic  finite  elements 
were  used  to  model  combustion  chamber  configurations 
(4,10,11).  These  studies  showed  good  agreement  between 
acoustic  finite  element  model  predictions  and  experimental 
results  for  various  shapes  of  combustion  chambers. 

Hickling  et  al  (4)  proposed  that  additional 
investigation  should  be  done  on  the  nodal  patterns  of  the 
cavity  resonances  within  the  cylinder  during  a  cycle,  or,  in 
other  words,  as  a  function  of  cylinder  depth.  Craggs  et  al 
(11)  investigated  the  cavity  resonances  and  mode  shapes  for 
a  Chevette  engine  combustion  chamber  as  a  function  of 
cylinder  depth.  This  provided  the  impetus  for  this  research. 

1.3  Purpose  of  the  Research 

The  purpose  of  this  research  is  to  examine  the  cavity 

resonances  of  three  axisymmetric  combustion  chamber-cylinder 

configurations,  see  Figure  1.2,  as  a  function  of  cylinder 

depth  using  acoustic  finite  elements.  The  combustion 

chambers  have  the  same  volume  1  and  the  cylinders  have  the 

same  diameter.  One  combustion  chamber,  Figure  1.2a  ,  has  a 

larger  diameter  than  the  cylinder.  In  Figure  1.2b  the 

combustion  chamber  has  the  same  diameter  as  the  cylinder  and 

’The  volumes  of  the  combustion  chambers  were  selected  to 
depict  a  compression  ratio  =  9. 
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^  Comb.  Chamber 
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c) 


Comb.  Chamber 


Cylinder 


Figure  1.2  Combustion  Chamber-Cylinder  Configurations 
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in  Figure  1.2c  the  combustion  chamber  has  a  smaller  diameter 
than  the  cylinder.  In  each  of  the  cases  considered,  several 
of  the  natural  frequencies  for  different  cylinder  depths  are 
examined.  The  mode  shape  for  the  lowest  mode  is  examined  in 
detail  for  each  of  the  three  cases  at  a  small  cylinder 
depth.  The  form  of  the  mode  shape  is  needed  in  order  to 
obtain  the  best  location  for  knock  detectors.  This  research 
is  concentrated  on  the  lowest  mode  because  engine  knock  has 
been  shown  to  be  related  to  the  lowest  acoustical  cavity 
resonance . 


2.  ACOUSTIC  THEORY 


2.1  Basic  Acoustic  Theory 

Some  of  the  fundamental  concepts  which  are  discussed 
throughout  this  report  are  those  associated  with  standing 
waves  of  sound  in  three-dimensional  enclosures.  An  excellent 
description  of  sound  in  enclosures  is  given  by  Morse  (5).  A 
brief  summary  of  some  of  the  acoustical  terminology  follows. 

Acoustic  pressure  is  defined  as  the  pressure  in  excess 
of  the  ambient,  see  Figure  2.1.  Zero  acoustic  pressure 
indicates  that  there  is  no  excess  pressure.  In  acoustic 
terms,  a  "nodal  surface"  is  a  three-dimensional  surface 
where  the  acoustic  pressure  is  zero. 

Sound  is  defined  as  fluctuations  in  the  acoustic 
pressure  due  to  some  vibrating  source.  The  vibrating  source 
creates  compressions  and  rarefactions  in  the  medium,  in  this 
case  air,  which  are  defined  as  wave  motion.  The  speed  of 
propagation  of  sound  waves,  in  a  given  medium,  is  defined  as 
the  velocity  of  sound.  It  is  a  function  of  the  density, 
pressure  and  temperature  of  the  medium.  This  study 
concentrates  on  sound  waves  in  air.  The  velocity  of  sound  in 
air  is  given  the  symbol  c: 


c  =  /  kRT 
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Figure  2.1  'Acoustic  Pressure 
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where 

k  =  the  ratio  of  specific  heats  of  air  =  1.4 

R  =  the  gas  constant  for  air  =  287  J/kg-°K 

T  =  the  absolute  temperature  =  °C  +  273  =  °K 

For  this  analysis,  the  results  are  given  at  T  =  25  °C. 

Therefore,  c  is  given  the  value  346  m/sec. 

Three-dimensional  enclosures  containing  air  have  normal 
modes  of  vibration.  If  a  sound  source  is  placed  within  a 
three-dimensional  enclosure,  it  will  excite  these  modes  of 
vibration  similar  to  the  excitation  of  standing  waves  on  a 
string.  These  sound  waves  are  a  function  of  space  (x,y,z) 
and  time  (t).  The  normal  modes  of  vibration  of  air  in  an 
enclosure  are  variations  in  the  acoustic  pressure  as  a 
function  of  space  and  time. 

For  this  analysis,  it  is  assumed  that  the  enclosures 
have  perfectly  smooth  rigid  walls.  The  boundary  condition 
for  a  rigid  wall  is  such  that  air  velocity  perpendicular  to 
the  wall  is  zero.  The  wave  equation  in  rectangular 
co-ordinates  (5)  is  given  as: 


2  2  2 

91P  +  +  OL 

2  2  2 
ax  9y  9z 


2 

l  ai 
?  at2 


(2.1.1) 


This  is  the  partial  differential  equation  of  the  acoustic 
pressure  p(x,y,z,t).  The  rigid  wall  boundary  condition  is 
given  as: 


(2.1.2) 
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For  the  analysis  of  the  combustion  chamber-cylinder 
configurations,  the  variation  in  the  acoustic  pressure  as  a 
function  of  space  only  is  used.  The  resulting  equation  is 
called  the  Helmholtz  equation.  The  Helmholtz  equation  can  be 
obtained  from  equation  (2.1.1)  using  a  separation  of 
variables  technique.  Therefore,  the  differential  equation 
used  in  this  analysis  is  given  in  terms  of  the  Laplacian 
operator  as: 

V2p  +  (^)2  P  -  0  (2.1.3) 


with  the  same  boundary  condition  as  given  by  equation 
(2.1.2),  or : 


n  •  Vp  =  0  (2.1.4) 

Using  the  above  differential  equation  (2.1.3)  and  the 
given  boundary  condition  (2.1.4),  the  Helmholtz  equation  can 
be  solved  analytically  for  certain  enclosures.  The  acoustic 
pressure  p  can  be  found  in  closed  form  solution  as  a 
function  of  the  spatial  co-ordinates.  Morse  (5)  deals 
specifically  with  the  solution  for  rectangular  and 
cylindrical  enclosures. 

For  a  rectangular  enclosure,  the  Helmholtz  equation  in 
Cartesian  co-ordinates  is  given  as: 
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The  solution  as  given  by  Morse  (5)  with  a  harmonic  time 
variation  exp(icjt)  is: 
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where  the  origin  of  the  (x,y,z)  co-ordinates  is  located  at 
the  centre  of  the  enclosure.  Equation  (2.1.6)  describes  the 
acoustic  pressure  variation  as  a  function  of  the  spatial 
co-ordinates  (x,y,z).  The  values  <jx  ,xjy  ,  and  u>2  represent  the 
natural  circular  frequencies  of  the  waves  in  each  of  the 
three  directions. 

For  a  cylindrical  enclosure,  the  Helmholtz  equation  in 
cylindrical  co-ordinates  is  given  as: 


+  1|p +  1  A  + a!a+  {u))2 

r  8r  r  ae2  az2  c 


p  =  0 


(2.1.7) 


The  solution  as  given  by  Morse  (5)  with  a  harmonic  time 
variation  exp(iwt)  is: 


p  =  c?s  (m6)  cos 
K  sin  x  ' 


m 


raj  rh 
r 


(2.1.8) 


where  the  origin  is  located  in  the  centre  of  one  of  the 
circular  ends  of  the  cylinder.  Equation  (2.1.8)  describes 
the  acoustic  pressure  variation  as  a  function  of  the  spatial 
co-ordinates  (r,0,z).  The  values  uz  and  u>v  represent  the 
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natural  circular  frequencies  for  waves  in  each  of  the  three 
directions  (cjr  encompassing  both  the  r  and  8  directions). 

The  quantity  Jm  represents  the  Bessel  function  which  is 
described  in  detail  by  Morse  (5).  Since  the  Bessel  function 
is  not  used  explicitly  in  this  analysis,  the  reader  is 
referred  to  Morse  (5)  for  further  explanation  if  so 
required . 

Equations  (2.1.6)  and  (2.1.8)  give  the  acoustic 
pressure  variations  for  two  different  enclosures.  These 
acoustic  pressure  variations  have  natural  frequencies  and 
mode  shapes  associated  with  them.  For  this  analysis,  the 
combustion  chamber-cylinder  configurations  are  examined  as 
rigid  enclosures  and  the  natural  frequencies  are  referred  to 
as  "cavity  resonances".  The  mode  shapes  are  important  in 
establishing  the  location  of  zero  acoustic  pressure  or 
"nodal  surfaces".  For  the  exact  solutions  the  mode  shapes 
are  given  by  the  sin/cos  function  for  the  rectangular 
solution  (2.1.6)  and  by  a  cos/sin  and  Bessel  function  for 
the  cylindrical  enclosure  (2.1.8). 

At  this  point,  the  problem  of  determining  cavity 
resonances  for  simple  rigid  enclosures  has  been  described  by 
a  differential  equation  and  boundary  condition  in  terms  of 
acoustic  pressure.  The  solution  yields  natural  circular 
frequencies  and  mode  shapes  obtained  in  closed  form  for 
rectangular  and  cylindrical  enclosures. 

Another  important  quantity  in  acoustic  terms  is  the 
energy  associated  with  the  acoustic  pressure  waves.  The 


* 


14 


energy  of  the  sound  wave  can  be  broken  down  into  potential 
energy  and  kinetic  energy  (5).  Using  the  derivation  of 
energy  intensity  for  a  plane  wave,  the  average  potential  and 
kinetic  energies  for  a  three-dimensional  wave  are  as 
follows : 


Potential  Energy  =  - 

2  pc 


2. 

p  dv 

- 

vol  ume 


(2.1.9) 


Kinetic  Energy  = 


1 

2pa)2 


(vp)2 

4 

vol  ume 


dv 


(2.1.10) 


These  equations  represent  the  average  energy  of  a 
simple  harmonic  acoustic  pressure  wave  of  frequency  f  = 
cj/2 7r .  These  energies  together  represent  the  total  average 
energy  of  sound  waves  in  enclosures  with  rigid  walls  as 
there  are  no  dissipative  agencies  associated  with  the 
pressure  waves  in  rigid  enclosures.  The  kinetic  and 
potential  energies  take  on  a  significant  meaning  in  the 
finite  element  analysis. 

2.2  Acoustic  Theory  as  an  Approximation  to  Complex 
Geometries 

Acoustic  theory  is  extremely  important  in  understanding 
the  concept  of  cavity  resonances  of  enclosures.  It  is 
particularly  important  to  be  able  to  compare  the  results  of 
numerical  models  to  some  exact  solution,  to  establish  bounds 
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for  the  numerical  results.  The  object  of  this  analysis  is  to 
model  complex  shaped  combustion  chamber-cylinder 
configurations  of  which  there  is  no  known  closed  form 
solution,  to  obtain  cavity  resonances  and  mode  shapes.  The 
exact  solution  for  rectangular  and  cylindrical  enclosures 
can  be  used  to  determine  the  convergence  of  the  numerical 
method  by  modelling  these  enclosures  and  comparing  the 
results  to  the  exact  solution.  This  procedure  gives  an 
indication  as  to  the  accuracy  of  the  acoustic  finite  element 
models  used. 

It  is  for  this  reason  that  the  solution  of  the 
Helmholtz  equation  for  a  cylindrical  enclosure  is  used 
repeatedly  throughout  this  analysis,  even  for  the  combustion 
chambers  whose  radii  differ  from  that  of  the  cylinder.  The 
exact  solution  is  used  to  establish  upper  and  lower  bounds 
for  the  numerical  results.  In  the  three  cases  considered, 
the  combustion  chamber  and  cylinder  are  axisymmetric 
cylindrical  enclosures  joined  together  to  form  a  more 
complicated  geometry  (see  Figure  1.2).  At  top  dead  centre, 
when  the  cylinder  depth  is  negligible,  the  solution  should 
be  governed  by  the  radius  of  the  combustion  chamber.  This 
establishes  a  bound  for  the  numerical  results  based  on  the 
analytical  solution  for  a  cylindrical  enclosure  with  a 
radius  equal  to  that  of  the  combustion  chamber.  As  the 
cylinder  depth  increases,  such  that  the  cylinder  volume 
dominates,  the  solution  should  be  governed  by  the  radius  of 
the  cylinder.  This  establishes  the  second  bound  for  the 
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numerical  results.  In  the  first  case,  the  mode  shapes  for 
the  combustion  chamber  should  be  similar  to  those  for  a 
cylindrical  enclosure.  In  the  second  case,  because  the 
combustion  chamber  is  present  the  mode  shapes  for  a 
cylindrical  enclosure  cannot  be  used  as  an  approximation. 

The  effect  of  the  combustion  chamber  on  mode  shapes  of  the 
combined  geometry  is  unknown. 

Rectangular  and  cylindrical  enclosures  will  be  dealt 
with  in  some  detail  in  the  next  two  sections. 

2.3  On  Rectangular  Enclosures 

Using  equation  (2.1.5)  with  the  rigid  wall  boundary 
condition,  the  solution  for  the  acoustic  pressure  variation 
in  a  rectangular  enclosure  as  given  by  Morse  (5)  is  given  by 
equation  (2.1.6).  This  leads  to  the  following  characteristic 
values  of  the  natural  frequencies  in  Hertz  which  represent 
the  cavity  resonances: 
f x  =  n  Xc/2LX 

f  y  =  nyc/2Ly  (2.3.1) 

fz  =  nzc/2Lz  2 

This  solution  gives  the  characteristic  values  of  a 
rectangular  enclosure  with  spatial  orientation  and 
dimensions  as  shown  in  Figure  2.2.  The  lengths  Lx,Ly,  and  Lz 
used  to  calculate  the  natural  frequencies  are  as  indicated. 

2Note:  nx,  ny,  and  nz  =  0,  1,  2  .  .  .  When  n  is  even,  the 
cosine  function  for  the  acoustic  pressure  variation  is  used. 
When  n  is  odd,  the  sine  function  is  used  (5). 
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Figure  2.2  Dimensions  and  Co-ordinate  System  for  a 


Rectangular  Enclosure 
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The  frequency  of  a  given  mode  as  given  by  Morse  (5)  is: 

1/2 
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The  frequency 

f 

in  Hertz  can 

be 

] 


(2.3.2) 


natural  frequencies  u>x ,  o)y ,  and  oz  according  to  equation 
(2.3.2)  and  knowing  the  velocity  of  sound  in  air  at  a  given 
temperature.  A  characteristic  value  sometimes  used  to 
describe  cavity  resonances  independent  of  c  is  (co/c)2. 

The  mode  shapes  for  a  rectangular  enclosure  are  given 
by  Morse  (5).  A  detailed  description  of  the  modes  of  a 
rectangular  enclosure  is  not  given  here  as  they  are  not  used 
explicitly  in  this  analysis.  The  reader  is  referred  to  Morse 
(5)  for  a  discussion  of  axial,  oblique  and  tangential  waves 
which  are  combinations  of  the  aforementioned  cavity 
resonances.  The  natural  frequencies  of  rectangular 
enclosures  are  used  later  to  study  the  accuracy  of  the 
acoustic  finite  element  analysis  compared  to  acoustic 
theory . 


2.4  On  Cylindrical  Enclosures 

Using  equation  (2.1.7)  with  the  rigid  wall  boundary 
condition,  the  solution  for  the  acoustic  pressure  variation 
in  a  cylindrical  encosure  as  given  by  Morse  (5)  is  given  by 
equation  (2.1.8).  This  leads  to  the  following  characteristic 
values  of  the  frequency  in  Hertz: 
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f  z  =  n  2  c /2L z 

f  =  amnc/2a 


(2.4.1) 


This  solution  gives  the  natural  frequencies  of  a 
cylindrical  enclosure  with  spatial  orientation  and 
dimensions  as  shown  in  Figure  2.3.  The  length  Lz  and  radius 
’a’  given  in  the  equations  (2.4.1)  to  calculate  the 
frequencies  are  as  indicated.  The  values  of  amn  are 
solutions  to  the  equation 


(ttoi)  =  0 


(2.4.2) 


which  are  related  to  the  Bessel  function  described  in 
equation  (2.1.8).  The  procedure  by  which  the  values  of  amn 
are  obtained  is  explained  in  detail  by  Morse  (5).  A  table  of 
various  values  of  amn  can  be  found  in  the  same  reference. 

For  this  analysis,  amn  for  the  lowest  four  modes  is  used.  In 
Morse’s  derivation  of  amn  the  subscripts  ’m’  and  'n'  denote 
the  type  of  mode  considered:  ’m’  represents  the  number  of 
diametral  nodal  surfaces  and  ’n'  represents  the  number  of 
radial  nodal  surfaces.  A  similar  description  of  modes  is 
used  in  this  analysis  except  three  indices  are  used  such 
that  the  axial  modes  are  included. 

Before  proceeding  further,  a  description  of  the  various 
types  of  nodal  surfaces  is  required.  Figure  2.4  illustrates 
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the  three  basic  types  of  inodes  present  in  a  cylindrical 
enclosure.  Figure  2.4a  has  a  diametral  nodal  plane  which 
extends  throughout  the  cylinder  with  no  variation  in  the  z 
direction.  This  is  referred  to  as  a  circumferential  mode. 
Figure  2.4b  has  a  cylindrical  nodal  plane  which  extends 
throughout  the  cylinder  with  no  variation  in  the  z 
direction.  This  is  referred  to  as  a  radial  mode.  Figure  2.4c 
has  a  nodal  plane  in  the  centre  of  the  cylinder  in  the  z 
direction.  This  is  referred  to  as  an  axial  mode.  For  each  of 
these  modes,  the  acoustic  pressures  on  either  side  of  the 
nodal  surface  have  opposite  phases  similar  to  that  of 
standing  waves  in  a  closed  tube. 

The  three  types  of  modes  are  referred  to  throughout 
this  analysis  in  the  following  manner: 

Mode  (mnz) 


where 

m  =  number  of  diametral  nodal  surfaces  (circumferential 
mode ) 

n  =  number  of  radial  nodal  surfaces  (radial  mode) 
z  =  number  of  axial  nodal  surfaces  (axial  mode)3 

As  the  values  of  ’m' ,  'n*  and  ’z’  can  each  take  on 
integer  values  from  0  to  combinations  of  the  different 
basic  types  of  modes  can  exist.  However,  in  this  analysis, 

3  Nodal  surfaces  are  surfaces  of  zero  acoustic  pressure. 
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which  concentrates  on  the  lowest  modes,  these  combined  modes 
are  not  considered.  Mode  (000)  is  the  rigid  enclosure  mode 
which  has  a  frequency  of  0  Hz  and  is  not  considered  to  be 
the  lowest  mode  of  concern. 

The  circumferential  and  radial  modes,  as  seen  from 
equation  (2.4.1)  are  dependent  on  the  radius  of  the 
cylindrical  enclosure.  The  axial  modes  are  dependent  on  the 
axial  length.  For  this  analysis,  the  modes  examined  are  the 
first  three  circumferential  modes,  the  first  radial  mode  and 
the  first  axial  mode.  Figure  2.5  shows  cross  sectional  views 
of  these  modes  in  their  relative  order  of  magnitude.  The 
nodal  surface  in  three  dimensions  is  projected  as  a  dashed 
line  in  the  plane  cross  section. 

For  the  cylindrical  enclosure,  the  lowest  mode  is  the 
first  circumferential  mode  (100)  with  one  nodal  diameter  as 
shown  in  Figure  2.5a.  Figure  2.5a  shows  a  pair  of  modes  for 
this  frequency  which  are  orthogonal  to  each  other.  From  the 
analysis  given  by  Morse  (5)  the  circumferential  modes  occur 
in  orthogonal  pairs  at  the  same  frequency.  Physically,  this 
would  mean  that  the  location  of  the  diametral  nodal  line  is 
not  fixed. 

The  next  lowest  mode  is  the  second  circumferential  mode 
(200)  which  has  two  nodal  diameters.  Mode  (200)  and  its 
orthogonal  counterpart  are  illustrated  in  Figure  2.5b.  The 
third  mode  is  the  first  radial  mode  (010)  which  has  one 
radial  nodal  surface  as  shown  in  Figure  2.5c.  This  mode 
exists  without  an  orthogonal  counterpart.  Mode  (300)  is  the 


, 


■ 


24 


Mode 


Mode 


Mode 


Mode 


Figure  2.5  Cross  Sectional  Representation  of 


(100) 


(200) 


(OlO) 


(300) 


(OOl) 


the  First  Five 


Modes  in  a  Cylindrical  Enclosure 


25 


fourth  mode  which  is  a  circumferential  mode  with  three  nodal 
diameters.  This  mode  exists  as  a  pair  of  orthogonal  modes, 
of  the  same  frequency  similar  to  the  other  two 
circumferential  modes.  The  axial  mode  (001)  is  placed  fifth 
in  order  of  magnitude.  However,  this  mode  is  dependent  on 
the  axial  length  rather  than  the  radius.  Therefore,  its 
order  of  magnitude  varies  depending  on  the  length  of  the 
cylindrical  enclosure  considered.  This  becomes  important  in 
the  analysis  of  the  combustion  chamber-cylinder 
configurations  where  the  cylinder  depth  changes  such  that 
the  relative  order  of  the  modes  could  change.  Also,  in  more 
complicated  geometries  such  as  complex  combustion 
chamber-cylinder  configurations,  purely  radial, 
circumferential  or  axial  modes  may  no  longer  exist. 

At  this  point,  the  types  of  modes  which  exist  in  a 
cylindrical  enclosure  have  been  described.  It  is  also 
necessary  to  obtain  a  numerical  value  for  the  frequencies 
which  are  associated  with  these  modes.  Using  the  values  of 
amn  obtained  from  Morse  (5),  the  frequencies  in  Hertz  are 
given  as  follows  for  each  mode: 


Mode  (100)  f  =  0.5861  c/2a 


Mode  (200)  f  =  0.9722  c/2a 


(2.4.3) 


Mode  (010)  f  =  1.2197  c/2a 
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Mode  (300)  f  =  1.3373  c/2a 

Mode  (001)  f  =  c/2L 

These  values  are  provided  for  a  cylindrical  enclosure 
of  length  L  and  radius  'a'.  A  more  detailed  description  of 
modes  in  cylindrical  enclosures  is  given  by  Morse  (5).  For 
this  analysis  these  modes  or  some  combination  thereof  may  be 
expected  for  the  combustion  chamber-cylinder  configurations 


given . 


3.  FINITE  ELEMENT  THEORY 


3 . 1  Basics 

It  is  essential,  in  a  study  using  acoustic  finite 
elements,  to  have  some  knowledge  of  the  finite  element 
technique,  in  general.  There  are  basic  definitions  which  are 
used  throughout  the  literature  on  finite  elements  which 
apply  to  this  analysis.  Three  references  have  been  used  to 
obtain  information  on  the  finite  element  method,  in  general 
(12,13,14).  Much  of  the  work  in  these  references  is  credited 
to  various  other  authors. 

It  has  been  shown  earlier  that  a  continuum  such  as  air 
in  an  enclosure  can  be  represented  by  a  governing  equation 
and  boundary  condition  (equations  (2.1.3)  and  (2.1.4))  which 
often  cannot  be  solved  using  analytical  methods  due  to  some 
irregular  geometry.  The  idea  behind  finite  element  modelling 
of  a  continuum  is  to  approximate  the  solution  domain  of  the 
irregular  geometry  using  a  number  of  smaller  regular  shaped 
geometries  (subdomains)  such  as  cuboids,  tetrahedrons  or 
pyramids,  and  thereby  approximate  the  solution  to  the 
governing  equation  and  boundary  condition  using  these 
subdomains.  These  smaller  regular  shaped  entities  are  called 
elements . 

In  a  continuum  such  as  air  in  a  rigid  enclosure  the 
unknown  solution  is  the  acoustic  pressure  variation 
throughout  the  enclosure.  A  continuum  such  as  this  is 
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composed  of  a  infinite  number  of  unknowns.  Using  finite 
elements,  a  continuum  can  be  represented  by  an  assemblage  of 
elements  which  contain  a  finite  number  of  unknowns  (12). 

The  element  is  assigned  a  given  number  of  nodes  which 
represent  the  values  of  the  unknown  solution.  For  example, 
if  a  cuboid  element  is  selected  with  nodes  at  each  corner, 
eight  values  for  the  unknown  variable  can  be  found  if  one 
unknown  per  node  is  sought.  In  the  analysis  of  combustion 
chamber-cylinder  configurations,  a  knowledge  of  the  acoustic 
pressure  at  various  nodes  provides  an  approximate  solution 
to  the  problem.  In  finite  element  terminology,  in  this  case, 
the  acoustic  pressure  p,  which  is  a  scalar  quantity,  is  the 
field  variable. 

There  are  numerous  types  of  elements  presented  in  the 
literature  on  finite  elements.  Pafec  (15)  and  Zienkiewicz 
(12)  give  excellent  references  on  the  types  of  elements  that 
can  be  used  for  different  types  of  problems.  The  type  of 
element  which  is  selected  is  usually  dependent  on  the 
desired  accuracy  of  the  solution,  the  geometry  under 
consideration  and  even  the  availability  of  computer  storage 
in  solving  the  final  matrix  equations.  Often,  an  element  is 
used  from  the  literature  because  its  behaviour  has  been 
extensively  studied.  An  element  is  said  to  be 
’isoparametric'  if  it  is  capable  of  being  distorted  in  a 
particular  fashion.  If  the  shape  function  used  to  distort  an 
element  is  the  same  function  used  to  approximate  the  field 
variable,  the  element  is  isoparametric. 
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Once  the  element  type  is  chosen,  an  interpolation 
function  is  selected  to  satisfy  certain  continuity 
requirements.  The  interpolation  function  is  the  assumed 
approximation  of  the  field  variable  over  the  element 
expressed  in  terms  of  the  nodal  values  of  the  field 
variable.  In  this  analysis,  the  interpolation  function  is  a 
polynomial  which  is  easy  to  integrate  and  differentiate.  The 
degree  of  the  interpolation  polynomial  depends  on  the  number 
of  nodes  as  well  as  the  number  of  unknowns  per  node  to 
satisfy  continuity.  Zienkiewicz  (12)  and  Heubner  (13)  give 
detailed  information  from  earlier  work  by  other 
investigators  on  commonly  used  finite  elements  and  their 
interpolation  functions.  The  continuity  requirements  are 
described  in  these  references  on  finite  elements  for  the 
general  case.  Continuity  requirements  must  be  met  so  that 
the  approximate  solution  converges  to  the  exact  solution  as 
the  number  of  elements  increases. 

Before  discussing  continuity  requirements  further,  the 
relationship  between  the  governing  equation  and  the 
individual  finite  element  must  be  found.  Generally, 
throughout  the  literature  on  finite  elements,  there  are  four 
methods  which  are  used  to  establish  the  finite  element 
equations  or  matrix  equations  which  represent  the  governing 
equation:  1)  the  direct  approach,  2)  Galerkin’s  method,  3) 
the  variational  approach  and  4)  the  energy  method  (13).  The 
method  used  for  deriving  the  matrix  equations  depends  on  the 
nature  of  the  problem.  The  variational  procedure  is 
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discussed  here  in  some  detail  to  give  a  brief  overview  of 
the  mathematical  nature  of  the  finite  element  solution. 

The  variational  procedure  for  deriving  the  finite 
element  equations  from  the  governing  equation  yields  a 
volume  integral  which  is  called  a  functional.  This 
functional  represents  the  governing  differential  equation 
and  boundary  condition.  In  terms  of  the  calculus  of 
variations,  if  this  functional  is  forced  to  have  a 
stationary  value,  the  unknown  function  which  yields  this 
stationary  value  is  the  approximate  solution  to  the 
differential  equation  and  boundary  condition. 

The  functional  is  made  stationary  by  taking  the  first 
variation  and  setting  this  equal  to  zero.  This  yields  a  set 
of  matrix  equations  for  the  solution  of  the  problem.  The 
matrix  equations  contain  volume  integrations  of  the  terms  of 
the  interpolation  polynomial  which  are  usually  evaluated 
using  numerical  integration  techniques. 

The  continuity  requirements  depend  on  the  order  of  the 
derivatives  appearing  in  the  functional.  The  following 
continuity  requirements  for  the  interpolation  function  for 

functionals  with  (r+1)th  derivatives  must  be  met: 

r 

1)  C  continuity  on  element  interfaces 
(r+l) 

2)  C  continuity  within  an  element 

In  addition  to  satisfying  the  continuity  requirements, 
the  numerical  integration  technique  must  be  selected  to 
evaluate  the  integrals  in  the  finite  element  equations 


* 


. 
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within  certain  limits  of  accuracy.  According  to  Cook  (14), 
on  the  basis  of  previous  work,  if  the  numerical  integration 
scheme  evaluates  the  volume  of  the  element  exactly,  this  is 
an  adequate  integration  scheme  to  assure  convergence  if  the 
continuity  requirements  have  been  satisfied. 

The  matrix  equations  for  a  finite  element  represent  the 
properties  of  an  element  according  to  the  governing  equation 
and  boundary  condition.  To  represent  the  entire  continuum, 
the  elements  are  assembled  together  knowing  that  the  field 
variable  is  the  same  for  common  nodes  of  elements  in  the 
assemblage . 

The  assembly  procedure  yields  a  set  of  matrix  equations 
which  can  be  solved  to  obtain  the  unknown  nodal  values  of 
the  field  variable.  This  is  the  basic  procedure  for 
obtaining  a  finite  element  solution  to  a  given  problem.  An 
excellent  reference  giving  detailed  steps  on  finite  element 
solutions  is  given  in  (13). 

3.2  Acoustic  Formulation  of  the  Finite  Element  Equations 

To  determine  the  cavity  resonances  of  combustion 
chamber-cylinder  configurations,  a  solution  to  the  Helmholtz 
equation  (2.1.3)  for  the  given  boundary  condition  (2.1.4) 
must  be  obtained.  The  analytical  solutions  for  rectangular 
and  cylindrical  enclosures  have  been  discussed  in  some 
detail.  A  technique  for  establishing  bounds  for  more  complex 
geometrical  configurations  has  been  provided  in  Section  2.2. 
However,  a  more  detailed  description  of  the  cavity 
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resonances  as  a  function  of  cylinder  depth  for  complex 
shaped  combustion  chamber-cylinder  configurations  has  not 
been  established. 

The  finite  element  method  has  been  described  as  a 
method  to  approximate  the  solution  to  certain  differential 
equations  and  boundary  conditions.  The  procedure  outlined  in 
Section  3.1  is  used  to  obtain  the  acoustic  finite  element 
equations. 

By  applying  variational  principles  to  equation  (2.1.3) 
and  satisfying  the  boundary  condition  equation  (2.1.4),  the 
following  functional  is  obtained  in  terms  of  the  acoustic 
pressure  p: 


i  (p) 


[ (vp)2  -  (£)2  p2]  dv 

vol  ume 


(3.2.1) 


This  functional  can  be  obtained  using  an  energy  approach  as 
indicated  by  Gladwell  (9). 

The  finite  element  equations  are  derived  from  the 
functional  for  one  element.  Often,  the  subscript  (e)  is  used 
to  denote  that  the  equations  represent  one  element;  however, 
this  notation  will  be  suspended  for  this  analysis.  The 
reader  should  keep  in  mind  that  the  resulting  matrix 
equations  represent  an  individual  acoustic  finite  element. 

In  order  to  find  a  stationary  value  of  the  functional 
(3.2.1),  an  approximate  solution  must  be  assumed  for  the 
element  (12).  This  approximate  solution  is  given  in  terms  of 
the  interpolation  polynomial  described  in  Section  3.1,  which 
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is  chosen  to  satisfy  the  continuity  requirements.  For  the 
acoustic  problem,  the  interpolation  function  or  the  assumed 
approximation  for  the  acoustic  pressure  variation  should 
have 

1)  C°  continuity  on  element  interfaces,  and 

2)  C1  continuity  throughout  the  element. 

C°  continuity  exists  if  the  field  variable,  acoustic 
pressure  p,  is  continuous  at  element  interfaces.  C1 
continuity  exists  if  the  first  derivatives  of  the  field 
variable  are  continous,  in  addition  to  the  field  variable 
itself.  This  is  the  criterion  which  is  used  to  select  the 
interpolation  polynomial  for  the  acoustic  finite  element 
formulation.  The  interpolation  polynomial  is  usually  treated 
generally  for  the  development  of  the  matrix  equations,  and 
is  referred  to  as  a  row  vector  LNJ .  The  approximation  for 
the  acoustic  pressure,  then,  takes  the  form: 

p  =  LNJ  W  (3.2.2) 

where  {p}  is  a  column  vector  representing  the  unknown 
acoustic  pressure  at  each  node  of  the  element. 

To  minimize  the  functional  I (p)  in  equation  (3.2.1), 
the  first  variation  of  the  functional  is  set  to  zero  as 
described  in  Section  3.1. 

3I(P)  -  o 

api 


(3.2.3) 


' 


E*i 


, 


34 


is  used  for  i  =  1,2,  .  .  .  m  where  'm'  is  the  number  of 
nodes  in  the  element. 

Substituting  equation  (3.2.2)  for  the  acoustic  pressure 
in  equation  (3.2.1)  and  performing  (3.2.3)  the  following 
matrix  equations  for  an  element  are  obtained: 

<[S]  -  (£)2  M)  {p>  ■  0 


where 


S 


i  j 


-3Ni  3TT  dH.  3TT 
Jx  3x  +  3y  3y 


vol  ume 


+ 


3N.  3N^- 
3z~  3z 


dv 


(N.N . )  dv 

vol  ume 


(3.2.5) 


These  are  the  finite  element  equations  for  an  acoustic 
element  in  Cartesian  co-ordinates  (x,y,z).  For  this 
analysis,  curved  elements  are  needed  to  represent  the  curved 
boundaries  of  the  combustion  chamber-cylinder 
configurations.  This  transformation  necessitates  formulating 
the  acoustic  finite  element  equations  in  terms  of  a  natural 
co-ordinate  system  (£,7?,$)'.  The  transformation  is  a  Jacobian 
transformation  of  a  curved  element  in  a  global  co-ordinate 
system  (x,y,z)  to  a  regular  shaped  element  in  a  natural 
co-ordinate  system  (£,7?,$)  (see  Figure  3.1)  .  This 
transformation  is  described  by  B.M.  Irons  in  (12).  To 
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Figure  3.1  Jacobian  Transformation 
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summarize,  the  global  co-ordinate  derivatives  are  found  in 
terms  of  the  natural  co-ordinates  in  the  following  manner: 


(3.2.6) 


The  volume  integral  is  also  transformed  in  the  following 
manner : 


dxdydz  =  det  [J]  d£dnd£  (3.2.7) 

In  the  above  analysis  [J]  is  the  Jacobian  matrix,  [J]'1 
is  the  inverse  of  the  Jacobian  matrix  and  det[j]  is  the 
determinant  of  the  Jacobian  matrix.  The  Jacobian  matrix  is 
formulated  in  the  following  manner: 


[J] 


3N. 

1 3Ni 

Z  3Ni 

as  xi 

3N- 

an  xi 

9N  • 
Z  3n 

3N . 

Y.  1 

3n 

l 


9C 


(3.2.8) 
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When  the  isoparametric  formulation  is  used,  as 
discussed  in  Section  3.1,  LNJ  ,  which  is  the  shape  function 
in  (3.2.6)  and  (3.2.8),  is  the  same  function  used  as  the 
interpolation  function  for  the  acoustic  pressure.  The 
summations  in  (3.2.8)  take  place  on  i  where  i  =  1,  2,  .  .  . 
m  and  'm'  is  the  number  of  nodes.  4  Applying  this 
transformation  to  the  matrices  in  (3.2.5)  gives  the 
following  matrix  expressions  in  terms  of  the  natural 
co-ordinate  system  (£,7?,$)  (see  Figure  3.1): 


[S]  ■  f  df  JT  Jdet[j])d5dr,dt 


vol ume 


(3.2.9) 


[P]  =  (LNJ  LNJ  det  dSdndc 

* 

vol ume 


In  this  case,  the  interpolation  function  of  the 
acoustic  pressure  and  the  shape  function  are  expressed  in 
terms  of  the  natural  co-ordinate  system  (£,77,$).  Equations 
(3.2.9)  represent  the  matrices  which  describe  an  acoustic 
finite  element  in  terms  of  the  interpolation  function  LN  J. 
The  interpolation  function  can  be  broken  down  in  the 
following  form: 


4  x,  y,  and  z  represent  the  global  co-ordinate  values  of  the 
nodes . 
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LNJ  ■  LFJ  [G]'1 


(3.2.10) 


The  interpolation  function  LNJ  is  expressed  as  the  product 
of  a  row  vector  LFJ  and  an  inverse  matrix  [G]  of  the  natural 
co-ordinates  (£,r?,$)  of  each  of  the  nodes  expressed  in  terms 
of  the  polynomial  LFJ .  This  derivation  of  the  interpolation 
function  is  explained  in  more  detail  in  (13).  In  the 
derivation  of  the  equations  for  a  particular  element,  (see 
Section  4.1),  the  transformations  are  more  easily 
understood.  Substituting  the  expression  for  LNJ  given  in 
equation  (3.2.10)  into  the  matrix  equations  (3.2.9)  gives 
the  following  result: 


vol ume 


det[J])  d£dr|d£ 


(3.2.11) 


[p]  =  f  (Lg'1JTLfJTLfJ  [G]_1det  [J]  d?dnd? 


vol  ume 


These  matrices  represent  the  acoustic  finite  element 
properties  for  an  individual  acoustic  element.  Matrix  [S] 
represents  a  volume  integral  of  the  square  of  the  gradient 
of  the  acoustic  pressure.  Matrix  [P]  represents  a  volume 
integral  of  the  square  of  the  acoustic  pressure.  Going  back 
to  the  energy  considerations  of  a  sound  wave  in  Section  2.1, 
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the  [S]  matrix  represents  the  kinetic  energy  of  an  elemental 
enclosure  and  the  [P]  matrix  represents  the  potential  energy 
of  an  elemental  enclosure.  This  result  was  obtained  directly 
by  Gladwell  (9)  using  energy  considerations  only. 

Once  the  element  matrices  are  obtained,  the  elements 
can  be  assembled  together  summing  up  the  contributions  of 
kinetic  energy  and  potential  energy  of  common  nodes.  In  this 
analysis,  a  computer  subroutine  is  used  to  assemble  elements 
together.  Once  the  elements  are  assembled  together  a  kinetic 
energy  matrix  and  potential  energy  matrix  are  obtained  which 
represent  the  entire  enclosure  with  rigid  boundaries.  The 
matrix  equation  which  is  solved  to  obtain  the  acoustic 
cavity  resonances  (eigenvalues)  of  the  enclosure  and  their 
associated  mode  shapes  (eigenvectors)  is: 

([S]  -  (^)2  [P])  {p}  =  o  (3.2.12) 

This  is  an  eigenvalue  problem  which  is  solved  using  a 
standard  eigenvalue  subroutine.  The  eigenvalues  are  obtained 
in  terms  of  the  non-dimensionalized  natural  frequency 
(<j/c)2.  The  order  of  the  square  matrices  [S]  and  [P]  is 
dependent  on  the  number  of  nodes  describing  the  total 
enclosure.  This  is  the  number  of  degrees  of  freedom  of  the 
solution  which  determines  the  number  of  eigenvalues  and 
eigenvectors . 
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4.  DEVELOPMENT  OF  ACOUSTIC  ELEMENT  PR 15 


4.1  Finite  Element  Equations  for  PR15 

For  the  analysis  of  the  combustion  chamber-cylinder 
configurations,  a  three-dimensional  isoparametric  prism 
element  was  chosen  to  model  the  enclosures.  The  element  is 
shown  in  Figure  4.1  ,  being  referred  to  as  element  PR15. 

This  element  has  fifteen  nodes,  or  in  other  words,  fifteen 
degrees  of  freedom.  The  element  was  chosen  because  its 
isoparametric  form  can  easily  model  a  pie  shape.5 

The  interpolation  function  for  element  PR15  was  derived 
similar  to  the  procedure  given  by  equation  (3.2.10).  The 
interpolation  function  was  not  used  directly  in  computing 
the  finite  element  equations,  but  rather  the  technique 
described  in  Section  3.2  was  used  where 

|_NJ  =  |_F_[  [G]"1  (4.1.1) 

The  polynomial  *-FJ  was  one  used  by  Pafec  (15).  LFJ  is  a  row 
vector  expressed  in  terms  of  the  natural  co-ordinates 
(£,77,$)  as  follows: 

|_FJ  =  l_i  S  n  S2  Sn  n2  c  nc  Sc  Snc  c2  S2c  n2c  Sc2  nc2J  (4.1.2) 


5  This  is  a  desirable  characteristic  for  modelling 
three-dimensional  axisymmetric  cylinders. 
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Limits  of  £  ;  0  <  £5  1 

Limits  of  rj ;  0  <  17  <1 

Limits  of  £• ;  ”1  <  f  <  1 


Figure  4.1  Element  PR15  in  Natural  Co-ordinates 


42 


Figure  4.1  shows  element  PR15  in  the  (£,77,$) 
co-ordinate  system.  The  nodes  are  numbered  in  a  specific 
order  such  that  the  the  [G]  matrix  given  in  equation  (4.1.2) 
was  evaluated  in  the  following  manner: 


[G]  = 


1  ni  ci2  cini  ^i2  ci  Vi  5ici  cinici  ci2  *1%  ni%  ^i2  Vi2 


i  c2  ... 


n2c2 


1  C15  n15  ^15  C15n15  ’  '  ‘ 


n15C15 


(4.1.3) 


The  subscripts  of  the  natural  co-ordinates  in  the  [G]  matrix 
correspond  to  the  node  numbers  in  Figure  4.1.  This  system  of 
numbers  was  used  so  that  a  computer  subroutine  could  be 
written  to  evaluate  the  element  properties  for  numerous 
elements  by  repeatedly  using  this  scheme.  The  matrix  [G]  is 
a  15  x  15  matrix  containing  the  natural  co-ordinates.  The 
inverse  of  this  matrix  was  obtained  using  a  standard 
computer  subroutine. 

To  evaluate  the  element  matrices  given  in  equations 
(3.2.11),  which  consist  of  volume  integrals,  a  numerical 
integration  scheme  was  used.  Instead  of  the  traditional 
numerical  integration  schemes  for  a  triangle  presented  in 
references  on  finite  elements  (13),  a  modified  approach  was 
implemented. 
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The  integration  scheme  is  shown  in  Figure  4.2.  The 
element  triangular  face  was  subdivided  in  the  £r\  plane,  see 
Figure  4.2a,  into  16  triangles  of  equal  area.  The  value  of 
the  polynomial  LFJ  was  evaluated  at  the  centroid  of  each  of 
the  triangles  and  multiplied  by  a  weighting  factor.  The 
weighting  factor  was  equal  to  1/16  of  the  area  of  the 
triangular  face.  In  the  $  direction,  see  Figure  4.2b,  a 
3-point  Gaussian  quadrature  (16)  was  used  to  evaluate  the 
polynomial.  The  integration  scheme  allowed  the  volume 
integrals  for  the  finite  element  equations  to  be  evaluated 
numerically  with  computer  subroutines. 

The  values  of  the  polynomial  LFJ  were  calculated  at 
each  of  the  integration  points  on  three  faces  as  shown  in 
Figure  4.2b.  The  Jacobian  matrix  [J],  the  inverse  Jacobian 
and  the  determinant  of  the  Jacobian  were  also  evaluated  at 
each  of  the  integration  points.  Standard  computer 
subroutines  are  used  to  evaluate  the  inverse  matrix  and 
determinant.  The  values  of  the  parameters  in  the  volume 
integral  of  the  matrices  given  in  equation  (3.2.11)  were 
evaluated  at  each  integration  point  and  the  results  summed 
together.  The  volume  integrals  can  then  be  interpreted  as  a 
summation  in  the  following  form: 

[S]  -  {  ¥  [G-’]Tl|  (5,  n1  Ck)  f  U,  ni  ck>  H  n,  ck)JT 

.-i 


k=i  i=i 

ni  ^  ni  ck^  L^i  ni  ck^3n  ni  ck^ 


f  ^  Ck)J  [G]"  det  [J(^.  fli  Ck)l  WiWk 


16 


[P]  =  E  E  [G”1  ]T  |_F(5.  ^  ck)  J'  LF^i  hi  Ck)J  [G] 
k=1  1=1  det[J(Ci  T).  Ck)]W.Wk 


-1 


(4.1.4) 
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Figure  4.2  Numerical  Integration  Scheme  for  PR15 
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These  matrices  give  the  acoustic  properties  of  element 
PR15.  By  using  element  PR15  evaluated  in  natural 
co-ordinates  transformed  from  global  co-ordinates  (x,y,z), 
computer  subroutines  were  written  to  obtain  the  matrices 
given  in  equations  (4.1.4). 

4.2  Behaviour  of  Element  PR15 

An  acoustic  problem  is  set  up  in  a  global  co-ordinate 
system  (x,y,z).  To  illustrate  the  use  of  element  PR15,  a 
rectangular  enclosure  is  presented  composed  of  two  PR15 
elements.  Figure  4.3  shows  a  cuboid  with  sides  of  length  1 
mm. 

The  following  information  is  required  to  calculate  the 
cavity  resonances  and  mode  shapes  by  formulating  the  [S]  and 
[P]  matrices  according  to  equations  (4.1.4)  through  the  use 
of  the  subroutines  in  Appendix  1. 

1. (x,y,z)  co-ordinates  of  each  node,  according  to  the  number 
assigned  to  the  node  in  Figure  4.1. 

This  is  used  in  formulating  the  Jacobian  matrix  according  to 
equation  (3.2.8)  . 

2.  A  global  assembly  matrix  specifying  the  location  of  each 
of  the  nodes  for  the  assembly  process. 

This  is  done  according  to  the  numbering  system  prescribed  in 
the  development  of  PR15  in  Figure  4.1.  This  numbering  system 
is  used  because  the  computer  subroutines  evaluate  the 
element  in  the  natural  co-ordinate  (£,t?,$)  system  according 


to  this  scheme. 


•  ’ 
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Figure  4.3  Sample  Problem  of  a  Rectangular  Enclosure 
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The  computer  subroutines  were  written  to  numerically 
evaluate  the  [S]  and  [P]  matrices  for  an  element  through  two 
callable  subroutines.  The  [S]  and  [P]  matrices  for  the 
entire  enclosure  are  assembled  from  the  element  matrices 
using  the  global  assembly  matrix  in  an  assembly  subroutine. 
For  the  problem  of  the  rectangular  enclosure  the  following 
data  is  used: 


X 

Y 

Z 

1 

o 

• 

o 

o 

• 

o 

o 

• 

o 

2 

1.0 

o 

• 

o 

o 

• 
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3 
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• 
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1.0 

o 

• 

o 

4 

0.5 

o 

• 

o 

o 

• 

o 

5 

0.5 

in 

• 

o 

o 

• 

o 

6 

o 

• 

o 

0.5 

o 

• 

o 

7 

o 

• 

o 

o 

• 

o 

0.5 

8 

1.0 

o 

• 

o 

0.5 

9 

o 

• 

o 

1.0 

0.5 

10 

o 

• 

o 

o 

• 

o 

1.0 

1  1 

1.0 

o 

• 

o 

1.0 

12 

o 

• 

o 

1.0 

1.0 

13 

in 

• 

o 

o 

• 

o 

1.0 

14 

0.5 

o 

• 

cn 

1.0 

15 

o 

• 

o 

in 

• 

o 

1.0 

Global  Assembly  Matrix  (2,15) 
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14  20  16  17  18  15  10  12  1 1  173452 
22  1620  19  1821  13  11  12937658 

The  [S]  and  [P]  matrices  for  each  element  are  generated 
by  calling  the  subroutines  ES  and  EP.  The  subroutines  are 
included  in  Appendix  1. 

For  the  problem  under  consideration,  there  are  22 
degrees  of  freedom.  After  the  two  elements  have  been 
assembled,  the  problem  is  solved  by  finding  the  eigenvalue 
solution  to  equation  (3.2.12)  for  the  given  assembled 
matrices  [S]  and  [P],  The  are  22  eigenvalues  in  terms  of 
(cj/c)2  and  their  associated  eigenvectors.  The  results  can  be 
compared  to  acoustic  theory  by  examining  the  frequencies  for 
a  rectangular  enclosure  using  equations  (2.3.1).  A  rigid 
enclosure  mode  (cj/c)2  =  0.0000  Hz  should  appear  in  the 
solution.  The  results  for  the  rectangular  enclosure  are 
given  below  for  the  first  3  modes  compared  to  the  exact 
solution . 

Numerical  Exact 


(cj/c)2 


(cj/c  )  2 


00.0000 


0.00 


1 1.51343 


9.86960 


1 1.83771 


9.86960 
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11.99967  9.86960 

This  method  of  evaluating  the  natural  frequencies  of  an 
enclosure  can  be  used  to  examine  the  convergence  of  element 
PR15.  By  increasing  the  number  of  elements,  one  would  expect 
the  solution  to  converge  to  the  exact  acoustic  theory.  From 
earlier  considerations,  the  convergence  of  an  element  is 
dependent  on  the  integration  scheme,  in  addition  to  the 
continuity  requirements  of  the  interpolation  polynomial  (see 
Section  3.1).  The  volume  of  one  PR15  element  is  evaluated 
numerically  at  a  value  of  0.499997973.  This  is  slightly  less 
than  the  exact  value  of  0.500000000.  This  indicates  that  the 
integration  scheme  is  of  a  slightly  lower  order  and  the 
convergence  may  not  be  quite  as  expected,  ie.  the  solutions 
may  not  converge  monoton ically  from  above  (14). 

The  convergence  of  element  PR15  was  evaluated  by 
modelling  a  rectangular  enclosure  (similar  to  the  example) 
using  an  increasing  number  of  elements.  The  behaviour  of 
element  PR15  was  examined  by  increasing  the  number  of 
elements  first  in  the  z  direction  and  second  in  the  y 
direction.  Figures  4.4  &  4.5  show  the  configurations  used  to 
evaluate  the  convergence. 

For  the  z  direction,  Figure  4.4,  from  1  to  7  elements 
were  specified.  The  convergence  results  are  presented  in 
Figure  4.4  for  the  first  mode.  The  percentage  error  is  given 
in  terms  of  the  exact  solution  for  a  rectangular  enclosure 
of  the  given  dimensions.  The  solution  converges  from  above, 
but  actually  becomes  slightly  lower  than  the  exact  value  for 
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Figure  4.4 


Z  Convergence  of  PR15  for  lowest  mode 


51 


more  than  5  elements. 

For  the  y  direction,  Figure  4.5,  from  1  to  5  elements 
were  specified.  The  convergence  results  are  presented  in 
Figure  4.5  for  the  first  mode.  The  percentage  error  is 
given,  similar  to  the  z  convergence,  in  terms  of  the  exact 
solution.  The  convergence  in  the  y  direction  is  slightly 
different  from  that  in  the  z  direction  because  of  the 
integration  scheme  employed.  The  x  direction  convergence  is 
not  evaluated.  From  considerations  of  the  integration  scheme 
on  the  ?  face,  the  convergence  in  the  x  direction  should  be 
similar  to  that  in  the  y  direction. 
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Figure  4.5  Y  Convergence  of  PR15  for  lowest  mode 


5.  COMBUSTION  CHAMBER-CYLINDER  ANALYSIS 


5.1  Description  of  the  Problem 

The  three  combustion  chamber-cylinder  configurations 
considered  in  this  analysis  are  illustrated  in  Figure  5.1. 
All  three  configurations  have  a  cylinder  with  a  radius  of  50 
mm.  This  radius  was  arbitrarily  selected  to  be 
characteristic  of  a  spark  ignition  engine.  Combustion 
chamber  1,  Figure  5.1a,  has  a  radius  larger  than  the 
cylinder.  Combustion  chamber  2,  Figure  5.1b,  has  a  radius 
equal  to  that  of  the  cylinder.  Combustion  chamber  3,  Figure 
5.1c,  has  a  radius  smaller  than  that  of  the  cylinder.  In 
each  case,  both  combustion  chamber  and  cylinder  are 
axisymmetric  cylinders.  The  volumes  of  the  three  combustion 
chambers  are  approximately  equal. 

The  three  configurations  were  selected  to  demonstrate 
trends  in  the  cavity  resonances  as  a  function  of  cylinder 
depth  and  type  of  combustion  chamber.  The  depth  of  the 
cylinder  was  varied  from  1  mm  -  50  mm.  The  smallest  depth 
represents  a  top  dead  center  configuration  similar  to  that 
which  occurs  as  the  piston  approaches  the  top  of  its  stroke 
in  an  engine. 

Combustion  chamber  2  is  of  particular  interest  because 
the  radius  of  the  cylinder  is  equal  to  that  of  the 
combustion  chamber.  Therefore,  this  configuration  is  simply 
a  cylindrical  enclosure.  The  analysis  given  in  Section  2.4 
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Figure  5.1  Combustion  Chamber-Cylinder  Dimensions 
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applies  to  this  configuration  and  the  analytical  results  can 
be  used  to  examine  the  validity  of  the  numerical  technique 
(finite  element  method)  which  is  used. 

5.2  Analytical  Approximation  for  Bounds 

The  analysis  for  cylindrical  enclosures  can  also  be 
used  to  establish  bounds  for  the  numerical  results.  Based  on 
the  dimensions  presented  in  Figure  5.1,  the  bounds  were 
determined  in  terms  of  frequency.  These  values  were 
calculated  from  equations  (2.4.3)  for  modes  (100),  (200), 
(010),  and  (300)  substituting  for  the  radius  'a',  the  radius 
of  the  combustion  chamber  or  cylinder  as  indicated  earlier. 

‘  The  two  bounds  were  joined  together  by  a  curve  to  give  an 
idea  of  the  behaviour  of  the  cavity  resonances  using  the 
approximation  to  compare  to  the  numerical  results. 

The  rigid  enclosure  mode  was  not  considered  and  the 
axial  modes  were  not  examined,  at  this  point,  because  the 
axial  modes  are  a  function  of  cylinder  depth.  This  technique 
establishes  bounds  for  complex  combustion  chamber-cylinder 
configurations  and  could  even  possibly  be  used  as  an 
approximation . 

5.3  Finite  Element  Models  of  Axisymmetric  Combustion 
Chamber-Cylinder  Configurations 

The  three  finite  element  models  for  the  axisymmetric 

combustion  chamber-cylinder  configurations  are  presented  in 

‘The  radius  of  the  combustion  chamber  gives  one  bound  and 
radius  of  the  cylinder  gives  another. 
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Figure  5.2.  PR15  acoustic  finite  elements  were  used  to  model 
one  half  of  the  combustion  chamber-cylinder  configurations 
as  indicated  in  Figure  5.2.  The  details  of  the  actual  models 
are  not  presented,  that  is,  the  global  co-ordinates  and 
assembly  matrices.  The  standard  procedures  described  in 
Section  4.2  were  used  to  solve  the  problem. 

The  finite  element  models  are  slightly  different  for 
each  configuration.  Figure  5.2a  shows  a  three-dimensional 
half  chamber  model  composed  of  24  elements.  The  bottom  6 
elements  have  a  variable  depth  to  model  different  cylinder 
depths.  This  model  represents  combustion  chamber  1.  Figure 
5.2b  shows  a  three-dimensional  half  chamber  model  composed 
of  12  elements.  The  bottom  6  elements  have  a  variable  depth 
to  represent  the  different  cylinder  depths.  The  second  model 
represents  combustion  chamber  2.  Figure  5.2c  shows  a 
three-dimensional  half  chamber  model  composed  of  24 
elements,  where  the  bottom  18  elements  have  a  variable 
depth.  This  configuration  represents  combustion  chamber  3. 

A  computer  program  was  used  to  evaluate  the  [S]  and  [P] 
matrices  for  each  combustion  chamber-cylinder  configuration 
by  assembling  the  individual  elements  shown  in  Figure  5.2. 
The  eigenvalue  problem  given  by  equation  (3.2.12)  was  solved 
for  each  system,  for  a  given  cylinder  depth.  The  program  was 
rerun  to  obtain  several  different  cylinder  depths  giving  the 
cavity  resonances  of  each  configuration  as  a  function  of 
cylinder  depth. 
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The  results  for  each  combustion  chamber-cylinder 
configuration  were  obtained  in  terms  of  the 
non-dimensionalized  natural  frequencies  (cj/c)2,  which  are 
eigenvalues.  These  values  were  converted  to  frequency  in 
Hertz  using  c  =  346  m/sec.  The  computer  subroutine  used 
sorted  the  eigenvalues  in  order  of  magnitude  from  smallest 
to  largest.  For  each  configuration,  the  lowest  value  which 
was  obtained  was  (cj/c)  2  =  0.0000  which  represents  the  rigid 
enclosure  mode.  These  values  were  not  used  as  a  cavity 
resonance.  Hence,  the  lowest  mode  considered  was  the  lowest 
value  next  to  the  rigid  enclosure  mode. 

Because  the  theory  proposed  by  Draper  ( 1 )  suggests  that 
knock  is  related  to  the  lowest  cavity  resonance,  excluding 
the  rigid  enclosure  mode,  only  the  first  few  cavity 
resonances  were  examined.  In  terms  of  frequency,  a  range  of 
1,000-5,000  Hz  was  selected  to  examine  the  cavity  resonances 
for  cylinder  depths  ranging  from  1  mm  -  50  mm. 

Another  major  concern  in  this  analysis  is  the  mode 
shapes  associated  with  the  cavity  resonances.  From  the 
finite  element  results,  eigenvalues  and  eigenvectors  were 
obtained.  The  eigenvalues  were  the  cavity  resonances  as 
already  discussed  and  the  eigenvectors  were  the  mode  shapes. 
The  mode  shapes  give  an  indication  of  the  acoustic  pressure 
variation  within  the  modelled  enclosure  for  a  given 
frequency.  Values  of  the  acoustic  pressure  were  obtained  for 
each  node  in  the  finite  element  model  from  the  eigenvectors. 
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The  mode  shapes  from  the  finite  element  results  for  the 
three  different  combustion  chamber-cylinder  configurations 
were  examined  at  one  specific  depth  of  4.5  mm  at  the  cross 
section  at  the  end  of  the  cylinder.  This  depth  represents  a 
small  cylinder  depth  such  as  that  which  might  occur  near  top 
dead  center.  The  finite  element  results  show  the  acoustic 
pressure  variation  along  the  diameter  of  the  configuration 
for  a  cross  section  taken  at  the  bottom  of  the  cylinder.  The 
depth  selected  for  the  cylinder  and  the  cross  section 
examined  were  based  upon  foreknowledge  of  the  experimental 
setup  to  allow  for  comparisons  between  the  finite  element 
method  and  measured  values. 

5.4  Experimental  Models  of  Axisymmetric  Combustion 
Chamber-Cylinder  Configurations 

A  question  of  considerable  importance  in  any 
mathematical  analysis,  be  it  numerical  or  analytical,  is  the 
validity  of  the  analysis  in  comparison  to  experimental 
values.  It  is  for  this  reason  that  experimental  models  were 
constructed  to  simulate  the  combustion  chamber-cylinder 
geometries  at  various  cylinder  depths. 

The  experimental  models  are  schematically  shown  in 
Figure  5.3.  Photographs  of  the  models  follow  the  schematics 
in  Plates  1-  3.  The  experimental  models  consisted  of  200  mm 
x  200  mm  panels  of  plexiglass  with  holes  bored  through  the 
centers  of  inner  dimensions  comparable  to  the  radii  of  the 
combustion  chambers  (see  Figure  5.1)  and  cylinders.  The 
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Figure  5.3  Schematic  of  the  Experimental  Models 
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Plate  5.1  Experimental  Model  for  Combustion  Chamber  1 


Plate 


5.2  Experimental  Model 


for  Combustion  Chamber  2 


Plate  5.3  Experimental  Model  for  Combustion  Chamber  3 
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small  holes  located  mid-side  on  each  panel  represent  holes 
for  fastening  screws.  In  each  of  the  three  cases,  the  panels 
were  placed  together  in  different  combinations  to  obtain  the 
three  different  combustion  chambers  at  different  cylinder 
depths.  The  cylinder  portion  of  the  model  consisted  of 
several  plexiglass  panels  with  holes  of  radius  50  mm.  The 
panels  had  depths  of  =*4.5  mm,  =*9.0  mm  and  =*18.0  mm  which 
were  placed  together  to  obtain  cylinder  depths  ranging  from 
=*(4.5  mm  to  46.0  mm).  The  combustion  chamber  portion  of  the 
model  consisted  of  three  different  panels  of  radii  40  mm,  50 
mm  and  62.5  mm  with  depths  of  =*  9.9  mm,  =*  12.1  mm,  and  =* 

18.5  mm  respectively.  The  depth  of  the  combustion  chambers 
were  not  exactly  the  same  as  those  in  the  finite  element 
model  because  standard  thicknesses  of  plexiglass  material 
were  used.  Plexiglass  was  used  to  simulate  the  rigid 
boundary  condition.  The  panels  were  fastened  together  with 
screws  to  obtain  the  seven  cylinder  depths  for  the  three 
cases  considered. 

From  the  schematics  in  Figure  5.3,  it  can  be  seen  that 
the  inner  dimensions  of  the  plexiglass  models  represented 
the  combustion  chamber-cylinder  configurations  used  in  the 
finite  element  analysis  (see  Figure  5.1). 

An  200  mm  x  200  mm  plexiglass  panel7  with  a  25  mm  hole 
drilled  off  center  was  mounted  on  each  combustion  chamber  as 
shown  in  Figure  5.3.  This  allowed  a  25  mm  microphone  to  be 
placed  such  that  the  microphone  was  flush  mounted  in  the 
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combustion  chamber  slightly  off  center.  This  panel  was 
screwed  together  (fastening  holes  as  indicated  in  Figure 
5.4)  with  the  combustion  chamber  and  cylinder  plexiglass 
panels  (see  Figure  5.3)  to  obtain  various  cylinder  depths.  A 
300  mm  x  200  mm  plexiglass8  panel  as  shown  in  Figure  5.4 
with  a  6  mm  hole  drilled  in  the  center  was  placed  over  the 
end  of  the  cylinder.  This  panel  was  not  fastened  to  the 
combustion  chamber-cylinder  configuration  so  that  movement 
of  the  panel  was  possible.  A  6  mm  microphone  was  placed  in 
the  hole  such  that  the  microphone  was  flush  mounted  in  the 
cylinder.  The  source  and  receiver  panels  were  marked  TOP  and 
BOTTOM  to  indicate  the  orientation  used  during  the 
experiments.  The  models  were  turned  upside  down  to  conduct 
the  experiments  such  that  the  cylinder  was  above  the 
combustion  chamber.  This  was  done  purely  to  accomodate 
movement  of  the  receiving  panel. 

A  general  schematic  is  given  in  Figure  5.5  for  the 
experimental  models  of  the  combustion  chamber-cylinder 
configurations.  Each  of  the  three  configurations  can  be 
described  in  this  general  manner  showing  an  enclosure  of 
air.  The  experimental  models  were  set  up  such  that  there  was 
an  axisymmetric  cylindrical  plexiglass  enclosure 
representing  the  combustion  chamber  with  a  25  mm  microphone 
mounted  flush  in  the  chamber.  There  was  another  axisymmetric 
cylindrical  plexiglass  enclosure  joined  and  open  to  the 
combustion  chamber  to  represent  the  cylinder,  with  a  6  mm 
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Figure  5.4  Source  and  Receiving  Panels 
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microphone  mounted  flush  with  the  cylinder.  The  6  mm 
microphone  acted  as  a  receiver.  The  25  mm  microphone  acted 
as  a  transmitting  source.  By  exciting  the  air  inside  the 
combined  combustion  chamber-cylinder  configuration  with  the 
25  mm  source  and  measuring  the  response  with  the  6  mm 
receiver,  the  cavity  resonances  were  experimentally 
determined.  In  addition,  by  moving  the  location  of  the  6  mm 
receiver  on  the  surface  of  the  cylinder,  the  mode  shape  at 
that  cross  section  for  a  given  depth  of  4.5  mm  was 
determined . 

5.5  Experimental  Procedure 

The  evaluation  of  the  cavity  resonances  and  mode  shapes 
for  the  three  combustion  chamber-cylinder  configurations 
using  the  experimental  models  followed  the  same  basic 
procedure  for  each  configuration.  The  problem  was  examined 
in  terms  of  the  general  schematic  given  in  Figure  5.5.  There 
was  air  in  a  rigid  enclosure  excited  by  a  source.  The 
response  was  measured  with  a  receiver.  The  panels  were 
screwed  together  to  obtain  a  given  combustion 
chamber-cylinder  configuration  at  a  specific  depth. 

The  equipment  used  in  the  experimental  analysis  is 
shown  in  Figure  5.6.  The  source  was  a  25  mm  condenser 
microphone  excited  by  a  sinusoidal  sweep  of  frequencies 
using  a  function  generator.  The  25  mm  microphone  was  offset 
so  that  it  was  not  located  on  a  potential  nodal  surface.  The 
sweep  was  used  to  examine  the  response  of  the  enclosure  in 
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Figure  5.6  Experimental  Equipment 
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the  range  of  frequencies  from  1,000  -  5,000  Hz  (see  Section 
5.2).  The  response  was  measured  with  a  6  mm  condenser 
microphone.  The  response  from  the  6  mm  microphone  was  passed 
through  a  microphone  amplifier  via  a  preamplifier  into  a 
dual  channel  spectrum  analyzer  to  obtain  a  trace  of 
frequency  versus  amplitude  in  the  frequency  range. 
Approximately  1000  samples  were  analyzed.  These  traces  are 
called  frequency  spectra.  All  experiments  were  conducted  in 
an  anechoic  chamber  at  T  »  25  °C. 

The  frequency  spectra  were  obtained  for  each  of  three 
combustion  chamber-cylinder  configurations  at  seven  cylinder 
depths,  ranging  from  =(4.5  mm  -  46.0  mm).  Four  locations 
(see  Figure  5.7)  were  used  to  determine  the  frequency 
spectra  of  a  given  configuration  at  a  specific  depth.  The 
locations  are  shown  in  plan  view  relative  to  the  cylinder, 
because  the  receiving  panel  was  mounted  on  the  end  of  the 
cylinder  to  take  measurements.  These  four  locations 
(Configurations  I  -  IV)  were  chosen  such  that  no  modes  in 
the  given  frequency  range  were  missed  on  the  frequency 
spectra,  ie.  when  the  receiving  microphone  was  placed  on  a 
nodal  line.  A  nodal  line  exhibits  zero  acoustic  pressure  and 
therefore  no  reading  from  a  microphone  would  be  obtained.  By 
measuring  at  different  locations  and  examining  four 
frequency  spectra,  all  the  cavity  resonances,  indicated  by 
peaks  on  the  spectrum  trace,  in  the  given  frequency  range, 
should  have  been  present.  In  the  second  measurement  location 
(Configuration  II),  where  the  microphone  was  placed  in  the 
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center  of  the  cylinder  cross  section,  only  radial  or  axial 
modes  should  be  present.’ 

In  addition  to  assuring  that  all  cavity  resonances  were 
experimentally  determined,  the  four  measurement  locations 
gave  a  rough  idea  as  to  the  types  of  modes  on  the  cylinder 
cross  section  the  given  frequencies  represented,  by  relating 
the  position  of  the  microphone  to  the  mode  type. 

Figure  5.8  gives  sample  frequency  spectra  or  traces  for 
combustion  chamber  2  at  a  depth  of  4.5  mm.  These  four 
frequency  spectra  represent  measurements  conducted  at  the 
four  different  measurement  configurations  shown  in  Figure 
5.7.  Cavity  resonances  are  shown  as  peaks  on  the  spectrum 
traces.  Figure  5.8b  indicates  the  presence  of  one  mode  at  a 
frequency  of  4150.00  Hz,  where  the  microphone  is  located  at 
the  center  of  the  cylinder  cross  section.  The  other  modes 
are  likely  circumferential  modes  because  the  frequencies  do 
not  appear  in  configuration  II. 

The  frequency  traces  are  shown  for  combustion  chamber  2 
in  Figure  5.8  because  it  represents  a  simple  cylindrical 
enclosure.  The  behaviour  of  the  frequency  spectra  at 
different  locations  can  be  examined  in  regard  to  the  mode 
shapes  discussed  earlier  (see  Figure  2.5). 

At  a  depth  of  4.5  mm  based  on  the  analytical  solution, 

four  cavity  resonances,  or  peaks  on  the  spectrum  traces, 

should  be  apparent  in  the  traces.  The  following  modes  should 

be  experimentally  determined:  mode  (100)  f!  =  2028  Hz,  mode 

’These  modes  are  the  only  type  without  a  nodal  line  passing 
through  the  center  of  the  cross  section. 
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(200)  f2  =  3364  Hz,  mode  (010)  f3  =  4220  Hz  and  mode  (300) 
f 4  =  4626  Hz,  in  that  order.  Upon  examining  the  four 
frequency  spectra  taken  at  the  different  locations,  four 
cavity  resonances  are  indicated  at  approximately  the  same 
frequencies  given  from  the  analytical  analysis. 

Configuration  II  indicates  that  the  resonance  at 
approximately  4150  Hz  is  a  radial  mode  and  the  other  three 
are  circumferential  modes.  Circumferential  modes  have  a 
point  of  zero  acoustic  pressure  at  the  center  of  the 
cylinder  (see  Figure  2.5),  hence,  the  receiving  microphone 
at  configuration  II  did  not  pick  up  these  modes.10 

The  mode  shapes  for  the  lowest  cavity  resonances  for 
each  of  the  three  combustion  chamber-cylinder  configurations 
at  a  depth  of  4.5  mm  on  the  cross  section  at  the  end  of  the 
cylinder  were  obtained.  In  addition,  the  mode  shapes  were 
obtained  for  the  second  and  third  cavity  resonances  of 
combustion  chamber  2,  in  a  similar  manner,  to  compare  to 
analytical  theory  and  the  finite  element  results  because 
combustion  chamber  2  is  a  simple  cylindrical  enclosure. 

The  mode  shapes  were  obtained  experimentally  by 

measuring  the  acoustic  pressure  at  various  locations  on  the 

cross  section  at  the  end  of  the  cylinder  by  moving  the 

receiving  panel  to  various  locations.  A  grid  of  21  points 

was  used  as  shown  in  Figure  5.9.  The  microphone  was  moved  to 

each  of  the  21  points,  where  a  reading  of  the  pressure  was 

measured  in  volts  using  the  instantaneous  value  of  amplitude 

10Note:  There  is  some  variation  in  the  frequencies  for 
different  measurement  locations. 
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Figure  5.9  Experimental  Mode  Shape  Grid 
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on  the  spectrum  analyzer  for  each  point.  By  plotting  these 
voltages  at  various  locations,  the  mode  shape  was  obtained. 
Not  all  of  the  21  points  were  used  in  the  mode  shapes  given 
in  this  analysis  as  the  picture  became  too  complicated.  The 
pressure  variation  along  one  diameter  was  plotted  to 
determine  the  mode  shapes  and  the  nodal  line  of  the  cross 
section  was  indicated  as  a  dashed  line. 

A  given  mode  was  excited  by  setting  the  function 
generator  to  a  sine  wave  of  a  given  frequency  which  was 
obtained  from  the  experimental  cavity  resonance  analysis. 

The  acoustic  pressure  distribution  was  measured  for  a  given 
frquency  at  the  locations  indicated  in  Figure  5.9.  For, 
example,  there  was  a  cavity  resonance  at  f  =  2018.75  Hz  for 
combustion  chamber  2  at  a  depth  of  4.5  mm.  To  examine  the 
mode  shape  on  the  end  of  the  cylinder,  the  function 
generator  was  set  to  2031  Hz.11  The  acoustic  pressure  was 
measured  at  21  points.  The  results  were  plotted  on  the 
cylinder  cross  section  to  demonstrate  the  mode  shape.  The 
nodal  line  of  the  cross  section  was  drawn,  using  linear 
interpolation  where  necessary. 

By  repeating  this  procedure,  experimental  mode  shapes 

for  the  other  two  combustion  chamber-cylinder  configurations 

at  a  depth  of  4.5  mm  were  obtained.  Likewise,  the  mode 

shapes  for  the  second  and  third  cavity  resonances  of 

combustion  chamber  2  were  obtained  for  comparison  with 

analytical  theory.  Combustion  chamber  2  has  been  repeatedly 

1 1  This  value  was  slightly  higher  than  the  frequency  given 
in  the  cavity  resonance  analysis. 
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studied  in  more  detail  because  it  is  a  simple  cylindrical 
enclosure  and  analysis  of  these  results  gives  an  indication 
of  the  validity  of  the  mathematical  models  both  approximate 
(finite  element  results)  and  exact  (analytical  solution).  12 


,2Note:  This  method  does  not  guarantee  that  the  mode  shape 
is  the  maximum  amplitude  associated  with  the  cavity 
resonance . 
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6.  COMBUSTION  CHAMBER-CYLINDER  RESULTS 


6.1  Cavity  Resonances 

The  bounds  for  the  numerical  results  were  calculated  on 
the  basis  of  section  5.2  as  an  approximation.  Figures  6.1  - 
6.3  show  the  approximate  curves  for  the  three  configurations 
for  the  first  four  modes. 

The  acoustic  finite  element  results  for  combustion 

chamber  1  are  presented  in  Figure  6.4,  indicated  by  solid 

lines.  The  curves  are  presented  in  order  of  magnitude  such 

that  the  lowest  curve  is  the  lowest  cavity  resonance,  the 

second  curve  the  second  lowest  cavity  resonance  and  so  on. 

The  unusual  shape  of  the  curves  for  the  higher  resonances 

(2nd,  3rd  .  .  .)  for  large  cylinder  depths  suggests  that 

some  unusual  behaviour  is  occurring.  13  These  changes  in  the 

curves  are  due  to  the  presence  of  axial  modes  as  the  length 

of  the  cylinder  increases.  Therefore,  the  curves  presented 

in  order  of  magnitude  do  not  necessarily  represent  one 

specific  mode  except  for  the  lowest  cavity  resonance  which 

does  not  exhibit  this  unusual  behaviour.  The  results 

presented  for  the  higher  modes  are  not  as  accurate  as  the 

lowest  mode  partly  because  the  finite  element  results  lose 

accuracy  for  higher  modes  due  to  the  approximate  assumed 

,3It  should  be  noted  here,  that  the  curves  do  not  occur  in 
pairs  as  did  the  modes  for  a  cylindrical  enclosure.  This  is 
due  to  the  use  of  a  half  chamber  finite  element  model  which 
gives  only  one  of  any  orthogonal  pair.  This  does  not  affect 
this  analysis  as  the  pairs  of  frequencies  are  identical  for 
axisymmetric  cylinders  such  as  these. 
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Figure  6,4  Cavity  Resonances  for  Combustion  Chamber  1 
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solution  for  the  acoustic  pressure.  The  results  are 
presented  merely  to  give  an  indication  of  the  behaviour  of 
the  different  cavity  resonances  as  a  function  of  cylinder 
depth.  Accuracy  of  the  higher  modes  can  be  obtained  by  using 
more  elements  to  model  the  enclosure.  However,  the  computer 
solution  of  large  eigenvalue  problems  is  expensive.  Figure 
6.1,  which  shows  the  bounds  for  the  numerical  results  for 
combustion  chamber  1  agrees  well  with  Figure  6.4.  The 
frequencies  for  the  lowest  mode  increase  with  increasing 
depth  in  both  cases.  Figure  6.1  does  not  illustrate  the 
effect  of  the  axial  modes  as  only  the  radial  and 
circumferential  modes,  which  are  indicated  by  arrows,  were 
examined  to  determine  bounds.  From  the  standpoint  of  engine 
knock  at  small  cylinder  depths,  the  axial  modes  do  not 
appear  to  play  a  significant  role  for  this  configuration. 

The  finite  element  results  for  combustion  chamber  2  are 
shown  in  Figure  6.5,  indicated  by  solid  lines.  The  curves 
are  plotted  in  the  same  manner  as  those  for  combustion 
chamber  1,  in  order  of  magnitude.  These  results  are  of 
particular  interest  because  the  combustion  chamber-cylinder 
configuration  is  simply  an  axisymmetric  cylinder.  Therefore, 
the  analytical  solution  given  in  section  2.4  can  be  used. 

The  cavity  resonances  (exact  solution)  calculated  from 
equations  (2.3.3)  are  also  presented  in  Figure  6.5,  shown  by 
dashed  lines.  The  analytical  solution  gives  a  good 
indication  of  the  accuracy  of  the  finite  element  model. 
Excellent  agreement  was  obtained  between  the  finite  element 
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Figure  6.5  Cavity  Resonances  for  Combustion  Chamber  2 
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approximation  and  acoustic  theory  for  the  lowest  mode. 

The  analytical  curves  are  labeled  with  arrows 
indicating  the  mode  types  which  the  frequencies  represent 
based  on  the  analytical  solution.  Modes  (100),  (200),  (010), 
and  (300)  are  shown  in  that  order.  The  first  axial  mode  is 
also  shown  which  occurs  at  large  cylinder  depths.  The 
presence  of  the  axial  mode  indicates  why  the  second 
resonance  curve  changes  shape  such  that  mode  (200)  exists  as 
the  second  lowest  resonance  until  a  depth  of  approximately 
37  mm  and  then  the  second  lowest  resonance  is  an  axial  mode. 

The  analytical  solution  shows  that  the  finite  element 
solution  loses  accuracy  at  the  higher  modes  compared  to 
acoustic  theory.  The  second  resonance  curve  from  the  finite 
element  analysis  is  approximately  4.5%  higher  than  the 
analytical  solution  for  mode  (200)  until  a  depth  of  37  mm  is 
reached.  At  a  depth  of  37  mm  the  finite  element  solution  is 
very  close  to  acoustic  theory  for  mode  (001).  This  change  in 
accuracy  can  be  explained  by  examining  the  finite  element 
model,  Figure  5.2b,  more  closely.  There  are  six  elements 
used  in  the  circumferential  direction.  Hence,  based  on  the 
convergence  of  element  PR15,  finite  element  predictions  in 
the  circumferential  direction  should  be  extremely  good  for 
mode  (100),  and  reaonably  good  for  mode  (200)  with 
decreasing  accuracy  for  the  higher  circumferential  modes. 
There  is  one  element  in  the  radial  direction.  Hence,  the 
finite  element  predictions  for  first  radial  mode  (010)  will 
not  be  very  accurate.  There  are  two  elements  used  in  the 
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axial  direction.  Therefore,  finite  element  predictions  for 
the  first  axial  mde  (001)  are  quite  good.  Going  back  to 
Figure  6.5,  the  differences  between  analytical  theory  and 
the  finite  element  results  can  be  justified  by  the  finite 
element  model.  The  extensive  study  of  combustion  chamber  2 
is  very  useful  in  understanding  the  behaviour  of  the  finite 
element  models. 

Figure  6.5  can  also  be  compared  to  Figure  6.2  which 
shows  the  bounds  for  the  numerical  results.  However,  in  this 
case  the  bounds  are  merely  the  analytical  solution  to  the 
problem. 

The  finite  element  results  for  combustion  chamber  3  are 
shown  in  Figure  6.6,  indicated  by  solid  lines.  The  curves 
are  plotted,  similar  to  those  for  the  other  two 
configurations.  The  lowest  cavity  resonance  is  a  smooth 
curve.  The  higher  cavity  resonances  are  unusually  shaped,  a 
similar  behaviour  to  that  exhibited  by  the  higher  modes  for 
combustion  chambers  1  and  2.  This  phenomena  is  attributed  to 
the  axial  modes  which  have  a  dependence  on  cylinder  depth 
rather  than  the  radius. 

A  comparison  of  Figures  6.3  and  6.6,  shows  a  good 
agreement  between  the  bounds  for  the  numerical  results  and 
the  numerical  results  themselves  for  combustion  chamber  3. 
The  frequencies  for  the  different  modes  tend  to  be 
decreasing  with  increasing  cylinder  depth.  Figure  6.3  does 
not  illustrate  the  axial  modes.  Only  bounds  for  modes  (100), 
(200),  (010)  and  (300)  are  given. 
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The  lowest  cavity  resonance  in  each  of  three  cases  is  a 
continuous  smooth  curve.  This  seems  to  indicate  that  one 
type  of  mode  is  present.  As  has  been  previously  shown,  the 
main  focus  of  this  analysis  is  to  examine  the  lowest  cavity 
resonance  in  detail.  Figures  6.4  -  6.6  show  the  lowest 
cavity  resonances  for  three  types  of  axisymmetric  combustion 
chambers  of  the  same  volume  as  a  function  of  cylinder  depth. 

Spectrum  traces  were  experimentally  obtained  for  each 
combustion  chamber-cy 1 inder  configuration  at  seven  cylinder 
depths  ranging  from  M4.5  mm  -  46.0  mm)  to  compare  to  the 
numerical  results.  The  experimental  cavity  resonances  are 
plotted  in  Figures  6.4  -  6.6  along  with  the  finite  element 
results.  The  experimental  results  are  indicated  by  circles, 
such  that  the  radius  of  the  circle  shows  the  experimental 
variation  in  frequency  for  the  measurements.  The  frequency 
spectrum  analyzer,  a  dual  channel  FFT  device,  had  a 
resolution  of  6.25  Hz  on  the  range  used.  In  addition,  the 
temperature  inside  the  anechoic  chamber  varied  from  24.2  - 
26.2  °C  during  the  course  of  the  measurements.  Therefore, 
the  peaks  on  the  spectral  traces  varied  from  configuration  I 
to  IV  by  as  much  as  30  Hz.  ,4  Repeatability  tests  were 
conducted  on  several  configurations  to  substantiate  that  the 
spectral  peaks  were  valid. 

The  experimental  results  in  all  three  cases  agree  well 
with  the  finite  element  results  for  the  lowest  cavity 
resonance.  The  agreement  with  combustion  chamber  2  is 


1 4This  variation  is  not  justified  by  the  resolution  alone. 
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exceptionally  good  especially  compared  to  the  analytical 
solution . 

6.2  Mode  Shapes 

The  mode  shapes  are  presented  in  Figures  6.7  -  6.11  in 
pairs  to  compare  the  finite  element  mode  shapes  to  the 
experimental  mode  shapes.  15  Each  figure  shows  the 
comparison  for  a  given  combustion  chamber  for  a  given  mode 
at  a  depth  of  4.5  mm  on  the  cylinder  cross  section.  Figures 
6.7a  -  6.9a  show  the  lowest  modes  predicted  from  the  finite 
element  analysis  for  combustion  chambers  1,  2  and  3 
respectively.  The  finite  element  predictions  indicate  that 
these  modes  are  first  circumferential  modes  (100)  at  the 
depth  indicated.  An  acoustic  pressure  variation  across  the 
diameter  is  plotted  for  the  finite  element  results  and  the 
nodal  line  is  shown  as  a  dashed  line.  Only  the  pressure 
variation  along  one  diameter  was  plotted  because  the  grid 
points  for  the  models  did  not  coincide  with  the  experimental 
grid  points.  However,  a  distinct  comparison  can  be  made 
between  the  experimental  mode  shapes  and  the  finite  element 
mode  shapes  by  considering  the  dashed  line  in  each  case. 

The  validity  of  the  mode  shapes  can  be  determined  in  a 
similar  manner  to  that  done  for  the  cavity  resonances. 
Combustion  chamber  2  can  be  examined  in  comparison  to  the 
analytical  theory.  Again,  for  simplicity,  the  modes  were 

1  5The  experiments  were  conducted  such  that  the  combustion 
chamber-cylinder  models  were  turned  upside  down  with  the 
receiving  panel  on  top. 
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^Finite  Element  f=  1724  Hz 


- -  Nodal  Surface 


b)Experimental  f  =  1 72 5  Hz 


Figure  6.7  Mode  Shapes  for  Lowest  Cavity 


Resonance-Combustion  Chamber  1 
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a)Finite  Element  f= 2023  Hz 


- Nodal  Surface 


b)Experimental  f=  2031  Hz 


Figure  6.8  Mode  Shapes  for  Lowest  Cavity 
Resonance-Combustion  Chamber  2 
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a) Fi n ite  Element  f  =  2388  Hz 


Nodal  Surface 


b)  Experimental  f=  2375  Hz 


Figure  6.9  Mode  Shapes  for  Lowest  Cavity 


Resonance-Combustion  Chamber  3 
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a) F i n i t e  Element  f=  3513  Hz 


Nodal  Surface 


b)Experimental  f  =  336 9  Hz 


Figure  6.10  Mode  Shapes  for  Second  Cavity 

Resonance-Combustion  Chamber  2 
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a)Finite  Element  f  =  4724  Hz 


Nodal  Surface 


b)Experimental  f  =  4200  Hz 


Figure  6.11  Mode  Shapes  for  Third  Cavity 


Resonance-Combustion  Chamber  2 
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examined  at  depth  of  4.5  mm.  From  Figure  2.5,  there  should 
be  one  circumferential  mode  with  one  nodal  diameter.  Figure 
6.8  shows  both  experimentally  and  from  the  finite  element 
results  a  cross  sectional  view  with  one  nodal  diameter. 
Secondly,  there  should  be  a  circumferential  mode  with  two 
nodal  diameters.  Figure  6.10  shows  both  experimentally  and 
from  the  finite  element  results  a  cross  sectional  view  with 
two  nodal  diameters.  The  third  mode  should  be  a  radial  mode. 
Figure  6.11  shows  both  experimentally  and  from  the  finite 
element  results  a  cross  sectional  view  with  a  radial  nodal 
line.  In  the  case  of  the  radial  mode,  linear  interpolations 
were  performed  to  locate  the  zero  acoustic  pressure 
location.  Both  experimentally  and  with  the  finite  element 
analysis,  the  mode  shapes  were  symmetric  with  respect  to  the 
nodal  lines.  In  the  higher  modes  from  the  finite  element 
analysis,  mode  (010)  for  example,  the  mode  shapes  began  to 
lose  their  expected  shapes  compared  to  analytical  theory. 
This  is  not  surprising  as  the  cavity  resonances  were  not  as 
accurate  for  the  higher  modes.  Again,  because  of  the  concern 
with  the  lowest  mode,  the  mode  shapes  seem  to  be  extremely 
good . 

An  interesting  note  to  mention  regarding  the 
measurement  of  the  mode  shapes  is  in  regard  to  the  actual 
value  of  the  frequency  required  to  excite  a  given  mode.  As 
discussed  earlier,  the  frequency  of  excitation  was 
established  from  values  obtained  from  the  frequency 
spectrums  for  a  given  cavity  resonance.  During  the  course  of 
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the  experiments,  it  was  found  that  the  excitation  frequency 
could  be  varied  with  certain  limits  while  still  maintainiing 
one  specific  mode  shape.  Considering  combustion  chamber  2, 
for  example,  a  mode  shape  was  obtained  at  an  excitation  of  f 
=  2019  Hz  which  was  very  similar  to  that  obtained  at  2031 
Hz.  The  mode  shape  at  2019  Hz  was  not  included  because  there 
were  some  difficulties  with  the  function  generator 
maintaining  an  excitation  of  2019  Hz.  However,  it  is 
interesting  to  note  that  the  actual  frequency  associated 
with  a  given  mode  shape  experimentally  is  not  very  precise. 
In  other  words,  the  mode  shape  is  not  necessarily  the 
maximum  amplitude  mode  shape  associated  with  resonance.  This 
is  of  considerable  importance  in  some  of  the  modes  where 
there  was  considerable  variation  in  frequency  from  run  to 
run . 

From  the  analytical  analysis  of  standing  waves  in 
enclosures,  a  positive  amplitude  on  one  side  of  a  nodal 
surface  is  expected  and  a  negative  amplitude  on  the  other. 
The  negative  amplitude  indicates  a  phase  difference  which 
was  not  obtained  experimentally.  However,  the  phase 
difference  was  obtained  from  the  finite  element  analysis  as 
indicated  in  the  figures  . 
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7.  FURTHER  CONSIDERATIONS 

Axisymmetric  combustion  chambers  represent  only  a 
portion  of  the  complex  combustion  chambers  as  discussed 
earlier.  However,  the  finite  element  analysis  provided  in 
this  research  would  be  equally  valid  for  asymmetric 
combustion  chambers. 

As  an  extension  to  the  experimental  work  which  was  done 
on  the  axisymmetric  combustion  chambers,  one  asymmetric 
combustion  chamber  was  examined  at  a  depth  of  4.5  mm.  This 
combustion  chamber  was  similar  to  combustion  chamber  3 
(radius  =  40  mm)  except  that  the  combustion  chamber  was 
offset  as  shown  in  Figure  7.1. 

The  cavity  resonances  were  experimentally  obtained  in  a 
similar  manner  to  those  for  the  axisymmetric  combustion 
chambers.  The  spectrum  trace  for  configuration  I  is  shown  in 
Figure  7.2. 

The  results  for  the  lowest  cavity  resonance  indicate 
that  two  frequencies  exist  at  very  near  to  the  same  value. 
This  indicates  that  the  asymmetry  of  the  combustion  chamber 
separates  the  orthogonal  pairs  of  circumferential  modes 
(100)  into  two  separate  modes,  of  slightly  different 
f  requenc ies . 

The  experimental  mode  shapes  which  were  obtained  for 

the  two  lowest  cavity  resonances  are  presented  in  Figure 

7.3.  Figure  7.3b  at  2400  Hz  shows  a  circumferential  mode 

(100).  Figure  7.3a  at  2300  Hz  shows  a  combined  mode.  14 

1 ‘The  excitation  frequencies  were  chosen  slightly  higher 
than  and  slightly  lower  than  each  mode  to  uncouple  the 
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Figure  7.1  Asymmetric  Combustion  Chamber 
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Figure  7.2  Spectrum  Trace  for  Asymmetric  Combustion  Chamber 
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a)Experimental  f=  2300  Hz 


b)Experimental  f  =  2400  Hz 


1 4 (cont ' d) 

Figure  7.3  Mode  Shapes  for  Asymmetric  Combustion  Chamber 
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If  one  superimposes  two  (100)  modes  orthogonal  to  one 
another,  the  mode  shape  given  in  Figure  7.3a  is  similar  to 
the  combination  of  these  two  modes,  with  some  distortion  due 
to  the  asymmetry  of  the  combustion  chamber.  Evidently,  the 
experimental  mode  shape  in  Figure  7.3a  (representing  2312.50 
Hz)  was  not  uncoupled  from  the  other  mode  (representing 
2381 .25  Hz) . 

This  experiment  was  conducted  as  a  recommendation  for 
areas  of  further  research.  For  more  complex  asymmetrical 
combustion  chambers,  finite  element  analyses  could  be 
performed  and  compared  to  experimental  results  similar  to 
that  done  in  this  work.  The  mode  shapes  as  a  function  of 
cylinder  depth  could  be  studied  to  shed  more  light  on  cavity 
resonances  of  complex  shaped  enclosures. 


1 ‘ (cont ’ d )modes . 
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8.  CONCLUSIONS 


8.1  Cavity  Resonances 

The  cavity  resonance  results  for  the  three  combustion 
chamber-cylinder  configurations  as  a  function  of  cylinder 
depths  are  summarized  in  Figures  6.4  -  6.6.  For  each 
combustion  chamber,  there  was  good  agreement  between  the 
finite  element  predictions  and  the  experimental  results  for 
the  lowest  cavity  resonance.  The  difference  between  the 
finite  element  results  and  the  experimental  results 
increases  for  the  higher  modes.  This  is  confirmed 
part icular i ly  for  combustion  chamber  2  where  the  comparison 
is  made  to  acoustic  theory.  Acoustic  theory  agrees  very  well 
with  experimental  values  for  all  the  modes  considered  for 
combustion  chamber  2,  within  the  frequency  range.  The 
approximation  for  the  complex  geometries  using  the 
analytical  equations  for  the  first  four  modes  is  quite  good 
for  the  geometries  considered.  It  gives  a  reliable 
indication  of  the  behaviour  of  the  modes  as  the  cylinder 
depth  changes. 

8.2  Mode  Shapes 

The  experimental  mode  shapes  and  finite  element  mode 
shapes  are  summarized  in  Figures  6.7  -  6.11.  For  each 
combustion  chamber,  there  was  good  agreement  between  the 
finite  element  predictions  and  the  experimental  results  for 
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the  lowest  cavity  resonance  at  the  given  depth.  For  all 
three  combustion  chambers  at  the  cylinder  depth  considered, 
the  lowest  cavity  resonance  was  a  circumferential  mode  (100) 
with  one  nodal  diameter.  The  agreement  with  analytical 
theory  for  combustion  chamber  2  mode  shapes  validates  the 
usefulness  of  the  finite  element  method  for  predicting 
cavity  resonances  of  rigid  enclosures. 

8.3  For  Engine  Knock 

Although  this  analysis  considered  air  in  combustion 
chamber-cylinder  geometries,  the  results  could  be  quite 
useful  in  considering  engine  knock  if  one  knows  the 
properties  of  the  gas  in  the  actual  engine  cycle.  The 
observed  trends  in  the  cavity  resonances  as  a  function  of 
cylinder  depth  would  be  the  same  for  a  different  velocity  of 
sound  as  would  the  location  of  the  nodal  surfaces.  However, 
the  actual  values  of  the  frequencies  would  be  much  higher  in 
a  running  engine.  This  analysis  does  give  some  insight  into 
the  acoustic  nature  of  complex  shaped  enclosures  in  regard 
to  engine  knock  as  well  as  insight  into  the  natural 
frequencies  of  complex  shaped  rigid  enclosures,  in  general. 
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APPENDIX  1 


C  THIS  SUBROUTINE  CALCULATES  THE  P  MATRIX  FOR 
C  A  15  NODE  TRIANGULAR  ISOPARAMETRIC  ELEMENT 
C 

SUBROUTINE  PR 1 5EP ( EP , X , Y , Z ,M, N ) 

DIMENSION 

G( 15, 15) , GAXI (16) , GAET ( 16) ,GAZE(3) ,EP(M,M) ,WK ( 1 5 , 1 5) 
DIMENSION  FTF ( 15, 15) ,  GT  ( 15, 15) ,X( 15) ,Y( 15) ,Z( 15) fW1 (3) 
CALL  GINV1 5 (G ,M) 

CALL  NUMIN (GAXI , GAET , GAZE , N ,W , W 1 ) 

CALL  FTXF ( FTF , N , M , W , W 1 , GAXI , GAET , GAZE , G , X , Y , Z ) 

DO  1  I =  1  ,  M 
DO  1  J= 1 , M 
WK ( I , J ) =  0 . 0 
DO  1  K= 1 , M 

1  WK ( I , J ) =WK ( I , J ) +FTF (I , K ) *G ( K , J ) 

DO  2  1  =  1  ,M 

DO  2  J=  1  , M 

2  GT ( I , J ) =G ( J , I ) 

DO  3  1  =  1  ,M 

DO  3  J=  1  ,M 
EP ( I , J)=0.0 
DO  3  K= 1 , M 

3  EP ( I , J  )  =EP ( I , J ) +GT ( I , K ) *WK ( K  ,  J ) 

RETURN 

END 

THIS  SUBROUTINE  CALCULATES  THE  KERNEL  OF  THE 
P  MATRIX  FOR  A  15  NODE  TRIANGULAR  ISOPARAMETRIC 
ELEMENT 

SUBROUTI NE  FTXF ( FTF , N , M , W , W 1 , GAXI , GAET , GAZE , G , X , Y , Z ) 
REAL  FT( 15) ,FTF(M,M) ,GAXI (N) ,GAET(N) ,JAC(3,3) ,JACT(3,3) 
DIMENSION  X(M),Y(M),Z(M),G(M,M),W1(3) ,GAZE(3) 

DO  3  1=1 ,M 
DO  3  J= 1 ,M 
3  FTF ( I , J)=0.0 
DO  5  L=1 ,3 
DO  4  K= 1 , N 
FT( 1 )=1 .0 
FT ( 2 ) =GAXI (K) 

FT ( 3 ) =GAET ( K ) 

FT ( 4 ) =GAXI (K ) **2 
FT ( 5 ) =GAXI (K)*GAET(K) 

FT ( 6 ) =GAET ( K ) **2 
FT ( 7 ) =GAZE ( L ) 

FT ( 8 ) =GAET ( K ) *GAZE ( L ) 

FT ( 9 ) =GAXI ( K ) *GAZE ( L ) 

FT ( 1 0 ) =GAXI (K)*GAET(K)*GAZE(L) 

FT ( 1 1 )=GAZE(L)**2 
FT ( 12) =GAXI (K) **2*GAZE(L) 

FT( 1 3 ) =GAET ( K ) **2*GAZE ( L ) 

FT ( 14) =GAXI ( K ) *GAZE ( L ) **2 
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FT ( 15)=GAET(K)*GAZE(L)**2 

CALL  JACOB 3 ( JAC , JACT , DETJ , K , L , M , GAXI , GAET , GAZE , N , X , Y , Z , G ) 
DO  1  I  =  1  ,  M 
DO  1  J=1,M 

1  FTF ( I , J)=FTF(I , J)+FT(I ) *FT ( J ) *DETJ*W*W 1 (L) 

4  CONTINUE 
DO  2  I  =  1  ,  M 
DO  2  J= 1 ,M 

2  FTF (I , J)=FTF(I f J) 

5  CONTINUE 
RETURN 
END 

THIS  SUBROUTINE  SETS  THE  NATURAL  CO-ORDINATES 
FOR  THE  15  NODE  TRIANGULAR  ELEMENT 

SUBROUTINE  PR  1 5 ( XI , ET , ZE , M) 

DIMENSION  XI (M) , ET (M) , ZE (M) 

XI (  1  )=0.0 
XI (2  )  =  1 . 0 
XI ( 3 )  =  0 . 0 
XI ( 4 )  =  0 . 5 
XI  (  5  )  =  0 . 5 
XI (6)=0 . 0 
XI  ( 7  )  =  0 . 0 
XI  ( 8  )  =  1 .0 
XI  (9)  =  0.0 
XI  (  1  0 )  =  0 . 0 
XI ( 1 1 )  =  1 . 0 
XI ( 12)  =  0. 0 
XI  (  13)  =  0.5 
XI  (  1 4 ) =0 . 5 
XI  (  1  5  )  =  0 . 0 
ET ( 1 ) =0 . 0 
ET ( 2 ) =0 . 0 
ET ( 3 ) = 1 .0 
ET ( 4 ) =0 . 0 
ET ( 5 ) =0 . 5 
ET ( 6 ) =0 . 5 
ET ( 7 ) =0 . 0 
ET ( 8 ) =0 . 0 
ET ( 9 ) = 1 .0 
ET ( 1 0 )  =  0 . 0 
ET ( 1 1 ) =  0 . 0 
ET( 1 2 ) = 1 .0 
ET ( 1 3 ) =0 . 0 
ET ( 1 4 ) =0 . 5 
ET( 1 5 ) =0 . 5 
ZE ( 1 )=- 1 . 0 
ZE ( 2 ) =- 1 .0 
ZE ( 3 ) =- 1 .0 
ZE ( 4 ) =- 1 .0 
ZE ( 5 ) =- 1 .0 
ZE ( 6 ) =- 1 . 0 
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ZE ( 7 ) =0 . 0 
ZE ( 8 ) =0 . 0 
ZE ( 9 ) =0 . 0 
ZE ( 1 0 )  =  1 . 0 
ZE ( 1 1 )  =  1 . 0 
ZE ( 1 2 )  =  1 .0 
ZE ( 1 3 )  =  1 .0 
ZE ( 1 4  )  =  1 .0 
ZE ( 1 5 )  =  1 .0 
RETURN 
END 

THIS  SUBROUTINE  FORMS  THE  GINV  MATRIX  FOR  THE 
NATURAL  CO-ORDINATES  OF  THE  15  NODE  TRIANGULAR 
ELEMENT 

SUBROUTINE  GINV15(G,M) 

DIMENSION  G(M,M) ,  LK  ( 15) ,MK( 15) ,XI ( 15) , ET ( 15) , ZE ( 15) 
CALL  PR  1 5 (XI ,ETf  ZE,M) 

DO  1  I  =  1  ,  M 
G(I ,  1  )  =  1  .0 
G(I , 2)=XI (I ) 

G ( I , 3 ) =ET ( I ) 

G ( I ,4 )=XI (I ) **2 
G ( I , 5 ) =ET ( I ) *XI ( I ) 

G ( I f 6)=ET(I ) **2 
G ( I , 7 ) =ZE ( I ) 

G ( I , 8 ) =ET ( I ) *  ZE ( I ) 

G ( I , 9 ) =XI (I )*ZE(I ) 

G(I , 10)=XI (I ) *ET ( I ) *ZE ( I ) 

G ( I , 1 1 )  =  ZE ( I ) *  *  2 
G ( I f 1 2 ) =XI (I ) **2*ZE( I ) 

G ( I , 1 3 ) =ET ( I ) **2*ZE ( I ) 

G ( I , 1 4 ) =XI (I ) *  ZE ( I ) *  *  2 
1  G(I , 1 5 ) =ET ( I ) *  ZE ( I ) *  *  2 
CALL  MI NV ( G , M , D , LK , MK ) 

RETURN 

END 

THIS  SUBROUTINE  CALCULATES  THE  S  MATRIX  FOR 
A  15  NODE  TRIANGULAR  ISOPARAMETRIC  ELEMENT 

SUBROUTINE  PR  1 5ES ( ES , X , Y , Z , M , N ) 

DIMENSION 

G( 1 5  f 1 5) ,  GT( 15, 15) ,GAXI (16) ,GAET( 16) ,GAZE(3) ,ES(M,M) 
DIMENSION  DFTXDF ( 15,15),WK(15,15),W1(3) 

DIMENSION  X(M) ,Y(M) ,  Z  (M) 

CALL  G I N V 1 5 ( G , M ) 

CALL  NUMIN (GAXI ,GAET,GAZE,N,W,W1 ) 

CALL  DFTF( DFTXDF, N,M,W,W1 ,GAXI ,GAET,GAZE,G,X, Y,  Z) 

DO  1  1=1 , M 
DO  1  J= 1 , M 
WK ( I , J)=0.0 
DO  1  K= 1 , M 


. 
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1  WK ( I , J ) =WK ( I , J ) +DFTXDF ( I , K ) *G ( K , J ) 

DO  2  1  =  1  ,M 

DO  2  J= 1 , M 

2  GT (I , J ) =G ( J , I  ) 

DO  3  I = 1 , M 

DO  3  J= 1 , M 
ES ( I , J)=0.0 
DO  3  K= 1 ,M 

3  ES ( I ,  J)=ES(I , J)+GT(I , K ) *WK ( K , J ) 

RETURN 

END 

THIS  SUBROUTINE  CALCULATES  THE  KERNEL  OF  THE 
S  MATRIX  FOR  A  15  NODE  ISOPARAMETRIC  ELEMENT 

r 

SUBROUTINE  DFTF (DFTXDF ,N ,M,W,W1 ,GAXI , GAET , GAZE , G , X , Y , Z ) 
DIMENSION  DF ( 3 , 1 5 )  , WK (3,15), WK 1(3,15) , DFTXDF (M,M) 

REAL  GAXI (N) , GAET ( N ) ,G(M,M) , JAC ( 3 , 3 ) ,JACT(3,3) 
DIMENSION  X(15),Y(15),Z(15) , GAZE ( 3 ) ,W1 (3) 

DO  1  I = 1 , M 

DO  1  J= 1 , M 

1  DFTXDF (I ,  J )  =  0 . 0 

DO  5  L= 1 , 3 

DO  2  K= 1 , N 

DF( 1 , 1 )=0.0 

DF (  1 , 2  )  =  1 . 0 

DF (  1  , 3 ) =0 . 0 

DF ( 1 , 4)=2.0*GAXI (K) 

DF ( 1 , 5 ) =GAET ( K ) 

DF ( 1 , 6 ) =0 . 0 
DF ( 1 ,7)=0.0 
DF ( 1 , 8 ) =0 . 0 
DF ( 1 f  9)=GAZE(L) 

DF ( 1 , 1 0 ) =GAET ( K ) *GAZE ( L ) 

DF( 1 , 1  1  ) =0 . 0 

DF ( 1 , 1 2 ) =2 . 0*GAXI ( K ) *GAZE ( L ) 

DF ( 1 , 1 3 ) =0 . 0 

DF ( 1 , 1 4 ) =GAZE ( L ) **2 

DF ( 1 , 1 5 )  =  0 . 0 

DF ( 2 , 1 )  =  0 . 0 

DF(2f 2) =0.0 

DF ( 2 , 3 )  =  1 . 0 

DF ( 2 , 4  )  =  0 . 0 

DF (2,5) =GAXI (K) 

DF ( 2 , 6 ) =2 . 0*GAET ( K ) 

DF(2, 7) =0.0 
DF (2,8) =GAZE ( L ) 

DF (  2 , 9 ) =0 . 0 

DF ( 2 , 1 0 ) =GAXI (K)*GAZE(L) 

DF ( 2 , 1 1 ) =0 . 0 
DF (2, 12)=0.0 

DF ( 2 ,  1  3 ) =2 . 0*GAET ( K ) *GAZE ( L ) 

DF(2f 14) =0 . 0 

DF ( 2 , 15) =GAZE ( L ) **2 


. 
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DF ( 3 , 1 ) =0 . 0 
DF ( 3 , 2 ) =0 . 0 
DF (3,3)=0.0 
DF (3,4)=0.0 
DF ( 3 , 5 ) =0 . 0 
DF (  3 , 6  )  =  0 . 0 
DF ( 3  f  7  )  =  1 . 0 
DF (3,8) =GAET ( K ) 

DF (3,9) =GAXI (K) 

DF ( 3 , 1 0 ) =GAXI ( K ) *GAET ( K ) 

DF ( 3  ,  1  1  )  =2 . 0*GAZE ( L ) 

DF ( 3 , 1 2 ) =GAXI ( K ) **2 
DF ( 3 , 1 3 ) =GAET ( K ) **  2 
DF ( 3 , 1 4 ) =2 . 0*GAXI ( K ) *GAZE ( L ) 

DF ( 3 , 1 5) =2 . 0*GAET(K) *GAZE(L) 

CALL  JACOB 3 ( JAC , JACT , DETJ , K , L , M , GAXI , GAET , GAZE , N , X , Y , Z  ,  G ) 
DO  3  1=1,3 
DO  3  J= 1  , M 
WK (I , J  )  =  0 . 0 
DO  3  LL= 1 ,3 

3  WK ( I , J ) =WK ( I  ,  J ) + JAC ( I , LL ) *DF ( LL , J ) 

DO  6  1=1,3 

DO  6  J= 1  , M 
WK 1 (I , J)=0.0 
DO  6  LL= 1 , 3 

6  WK 1 ( I  ,  J  )  =WK  1  ( I  ,  J ) + JACT ( I , LL ) *WK ( LL , J ) 

DO  4  I  =  1  ,  M 

DO  4  J=  1  ,  M 
DO  4  LL= 1 , 3 

4  DFTXDF ( I , J ) =DFTXDF ( I , J ) +DF ( LL , I ) *WK 1 ( LL , J ) *DETJ *W*W 1 ( L ) 
2  CONTINUE 

DO  7  1=1 ,M 
DO  7  J= 1  , M 

7  DFTXDF (I ,J)=DFTXDF(I  ,J) 

5  CONTINUE 
RETURN 
END 


THIS  SUBROUTINE  FORMS  THE  JACOBIAN  FOR  A  3  DIMENSIONAL 
15  NODE  TRIANGULAR  ISOPARAMETRIC  ELEMENT 

SUBROUTINE 

JACOB3 ( JAC , JACT , DETJ , K , L , M , GAXI , GAET ,GAZE,N,X,Y,Z,G) 

REAL 

JAC (3, 3) , JACT (3,3) ,GAXI (N) , GAET (N ) , GAZE ( 3 ) , X (M) ,Y(M) ,Z(M) 
DIMENSION  WK1 ( 15) ,WK2( 15) ,WK3( 15) ,G(M,M) ,DFXI (15) , DFET ( 15) 
DIMENSION  DFZE ( 15) 

DIMENSION  LK ( 3 ) , MK ( 3 ) 

DO  4  1=1,3 
DO  4  J= 1 , 3 
4  JAC (I , J )  =  0 . 0 
DO  5  1=1  ,M 
WK 1 ( I ) =  0 . 0 
WK2 ( I ) =0 . 0 


5  WK3 ( I ) =0 . 0 

DFXI ( 1 )=0.0 

DFXI ( 2 ) = 1 .0 

DFXI (3)=0.0 

DFXI ( 4 ) =2 . 0*GAXI (K) 

DFXI ( 5 ) =GAET ( K ) 

DFXI ( 6 ) =0 . 0 
DFXI (7)=0.0 
DFXI ( 8 ) =0 . 0 
DFXI ( 9 ) =GAZE ( L ) 

DFXI ( 1 0 ) =GAET ( K ) *GAZE ( L ) 

DFXI ( 1 1 )=0.0 

DFXI ( 1 2 ) =2 . 0*GAXI ( K ) *GAZE ( L ) 

DFXI ( 1 3 ) =0 . 0 

DFXI ( 14 )=GAZE(L) **2 

DFXI ( 1 5 ) =0 . 0 

DFET ( 1 )=0.0 

DFET ( 2 ) =0 . 0 

DFET ( 3 ) = 1 .0 

DFET ( 4 ) =0 . 0 

DFET ( 5 ) =GAXI (K) 

DFET ( 6 ) =2 . 0*GAET ( K ) 

DFET ( 7 ) =0 . 0 
DFET ( 8 ) =GAZE ( L ) 

DFET ( 9 ) =0 . 0 

DFET (10) =GAXI ( K ) *GAZE ( L ) 

DFET ( 1 1 )=0.0 
DFET ( 1 2 ) =0 . 0 

DFET ( 1  3 ) =2 . 0*GAET ( K ) *GAZE ( L ) 

DFET ( 1 4 ) =0 . 0 

DFET ( 15)=GAZE(L)**2 

DFZE ( 1 )=0. 0 

DFZE ( 2 ) =0 . 0 

DFZE ( 3 ) =0 . 0 

DFZE ( 4 ) =0 . 0 

DFZE ( 5  )  =  0 . 0 

DFZE ( 6 ) =0 . 0 

DFZE ( 7 ) = 1 .0 

DFZE ( 8 ) =GAET ( K ) 

DFZE ( 9 ) =GAXI (K) 

DFZE ( 1 0 ) =GAXI ( K ) *GAET ( K ) 

DFZE ( 1 1 ) =2 . 0*GAZE ( L ) 

DFZE ( 1 2 ) =GAXI 00**2 

DFZE ( 1 3 ) =GAET ( K ) **2 

DFZE ( 1 4 ) =2 . 0*GAZE (L) *GAXI (K) 

DFZE ( 15)=2.0*GAZE(L)*GAET(K) 

DO  1  J= 1 ,M 
DO  1  1  =  1,  M 

WK1 ( J ) =WK 1 ( J ) +DFXI ( I ) *G ( I , J ) 
WK2 ( J ) =WK2 ( J ) +DFET ( I ) *G ( I , J) 

1  WK3 ( J ) =WK3 ( J ) +DFZE ( I ) *G ( I , J ) 
DO  2  1  =  1  ,M 

JAC (1,1) = JAC (1,1) +WK 1 ( I ) *X ( I ) 
JAC (2,1) = JAC (2,1) +WK2 ( I ) *X ( I ) 
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JAC (3,1) = JAC (3,1) +WK3 ( I ) *X ( I ) 

JAC (1,2) = JAC (1,2) +WK 1(1) *Y (I) 

JAC (2,2) = JAC (2,2) +WK2 ( I ) *Y ( I ) 

JAC (3,2) = JAC (3,2) +WK3 ( I ) *Y ( I ) 

JAC (1,3) = JAC (1,3) +WK 1 ( I ) *  Z ( I ) 

JAC (2,3) = JAC (2,3) +WK2 ( I ) *  Z ( I ) 

2  JAC (3,3) = JAC (3,3) +WK3 ( I ) *  Z ( I ) 

CALL  MI NV ( JAC , 3 , DET J , LK , MK ) 

DO  3  1=1,3 
DO  3  J=  1  , 3 

3  JACT ( I , J ) = JAC ( J , I ) 

RETURN 

END 

SUBROUTINE  NUMI N ( GAXI , GAET , GAZE , N , W , W 1 ) 

DIMENSION  GAXI (N) , GAET ( N ) ,GAZE(3) ,W1 (3) 

GAXI ( 1 )=1 .0/12.0 

GAXI ( 2 ) = 1 .0/6.0 

GAXI ( 3 ) = 1 .0/12.0 

GAXI ( 4 ) = 1 .0/3.0 

GAXI ( 5 ) = 1 .0/12.0 

GAXI ( 6 ) = 1 .0/6.0 

GAXI (7)=1 .0/3.0 

GAXI (8)=5. 0/12.0 

GAXI (9)=1 .0/12.0 

GAXI ( 1 0 ) = 1 .0/6.0 

GAXI ( 1 1 ) = 1 .0/3.0 

GAXI ( 1 2 )=5. 0/1 2 . 0 

GAXI (13) =7 . 0/1 2.0 

GAXI ( 1 4 ) =8 . 0/12.0 

GAXI ( 1 5 ) =7 . 0/12.0 

GAXI ( 1 6 ) = 1 0 . 0/12.0 

GAET ( 1 )= 10. 0/12.0 

GAET(2)=8. 0/12.0 

GAET(3)=7. 0/12.0 

GAET(4)=7. 0/12.0 

GAET(5)=5. 0/12.0 

GAET ( 6 ) = 1 .0/3.0 

GAET ( 7 ) = 1 .0/3.0 

GAET(8)=5. 0/12.0 

GAET ( 9 ) = 1 .0/6.0 

GAET ( 1 0 ) = 1 .0/12.0 

GAET ( 1 1 ) = 1 . 0/1 2.0 

GAET ( 12) =1 .0/6.0 

GAET ( 1 3 ) = 1 .0/3.0 

GAET ( 14)= 1 .0/6.0 

GAET ( 1 5 ) = 1 . 0/1 2.0 

GAET ( 16)= 1 .0/12.0 

W=1 .0/32.0 

GAZE ( 1 ) =0 . 0 

GAZE(2)=0. 7745966692 

GAZE ( 3 ) =~GAZE ( 2 ) 

W1 ( 1 )=8. 0/9.0 
W1 (2)=5. 0/9.0 
W1 (3)=W1 (2) 


. 


s.? 
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